diff --git a/src/basics.cc b/src/basics.cc index fc4ddd7..7bfe278 100644 --- a/src/basics.cc +++ b/src/basics.cc @@ -9,7 +9,10 @@ * Eryn Wells */ +#include #include +#include + #include "basics.h" #pragma mark - Vectors @@ -271,6 +274,221 @@ LinearCombination(const Double k1, const Vector3& v1, k1 * v1.z + k2 * v2.z + k3 * v3.z); } +#pragma mark - Matrices + +/* static */ Matrix4 +Matrix4::Zero() +{ + Matrix4 m; + memset(m.mCells, 0, 16 * sizeof(Double)); + return m; +} + + +/* static */ Matrix4 +Matrix4::Identity() +{ + Matrix4 m = Zero(); + for (int i = 0; i < 4; i++) { + m.mCells[i * 4 + i] = 1.0; + } + return m; +} + + +/* static */ Matrix4 +Matrix4::Translation(Double x, + Double y, + Double z) +{ + Matrix4 m = Identity(); + m.mCells[3] = x; + m.mCells[7] = y; + m.mCells[11] = z; + return m; +} + + +/* static */ Matrix4 +Matrix4::Rotation(Double x, + Double y, + Double z) +{ + Matrix4 m = Identity(); + + if (x == 0.0 && y == 0.0 && z == 0.0) { + /* No rotation, just return the identity matrix. */ + } else if (x != 0.0 && y == 0.0 && z == 0.0) { + /* + * Fill in m with values for an X rotation matrix. + * + * [1 0 0 0] + * [0 cos(x) -sin(x) 0] + * [0 sin(x) cos(x) 0] + * [0 0 0 1] + */ + Double cosX = std::cos(x); + Double sinX = std::sin(x); + m.mCells[5] = cosX; + m.mCells[6] = -sinX; + m.mCells[9] = sinX; + m.mCells[10] = cosX; + } else if (x == 0.0 && y != 0.0 && z == 0.0) { + /* + * Fill in m with values for a Y rotation matrix. + * + * [ cos(y) 0 sin(y) 0] + * [ 0 1 0 0] + * [-sin(y) 0 cos(y) 0] + * [ 0 0 0 1] + */ + Double cosY = std::cos(y); + Double sinY = std::sin(y); + m.mCells[0] = cosY; + m.mCells[2] = sinY; + m.mCells[8] = -sinY; + m.mCells[10] = cosY; + } else if (x == 0.0 && y == 0.0 && z != 0.0) { + /* + * Fill in m with values for a Z rotation matrix. + * + * [cos(z) -sin(z) 0 0] + * [sin(z) cos(z) 0 0] + * [ 0 0 1 0] + * [ 0 0 0 1] + */ + Double cosZ = std::cos(z); + Double sinZ = std::sin(z); + m.mCells[0] = cosZ; + m.mCells[1] = -sinZ; + m.mCells[4] = sinZ; + m.mCells[5] = cosZ; + } else { + /* + * TODO: Rotation in more than one dimension. So do a general rotation + * matrix. There's some magic way to do this with matrix multiplication + * that avoids gimbal lock. I should figure out how to do it properly. + */ + assert(0); + } + + return m; +} + + +/* + * Matrix4::Matrix4 -- + */ +Matrix4::Matrix4() + : mCells() +{ } + + +/* + * Matrix4::Matrix4 -- + */ +Matrix4::Matrix4(const Double cells[16]) + : mCells() +{ + memcpy(mCells, cells, 16 * sizeof(Double)); +} + + +/* + * Matrix4::Matrix4 -- + */ +Matrix4::Matrix4(const Matrix4& rhs) + : Matrix4(rhs.mCells) +{ } + + +/* + * Matrix4::operator() -- + */ +Double& +Matrix4::operator()(const unsigned int row, + const unsigned int col) +{ + assert(row < 4 && col < 4); + return mCells[4*row + col]; +} + + +/* + * Matrix4::operator* -- + */ +Matrix4 +Matrix4::operator*(const Double rhs) + const +{ + return Matrix4(*this) *= rhs; +} + + +/* + * Matrix4::operator*= -- + */ +Matrix4& +Matrix4::operator*=(const Double rhs) +{ + for (int i = 0; i < 16; i++) { + mCells[i] *= rhs; + } + return *this; +} + + +/* + * Matrix4::operator* -- + */ +Matrix4 +Matrix4::operator*(const Matrix4& rhs) + const +{ + return Matrix4(*this) *= rhs; +} + + +/* + * Matrix4::operator*= + */ +Matrix4& +Matrix4::operator*=(const Matrix4& rhs) +{ + Matrix4 lhs(*this); + for (int i = 0; i < 4; i++) { + for (int j = 0; j < 4; j++) { + /* Each cell is Sigma(k=0, 4)(lhs[ik] * rhs[kj]) */ + const int cell = i*4 + j; + mCells[cell] = 0.0; + for (int k = 0; k < 4; k++) { + mCells[cell] += lhs.mCells[i*4 + k] * rhs.mCells[k*4 + j]; + } + } + } + return *this; +} + + +/* + * Matrix4::CArray -- + */ +const Double* +Matrix4::CArray() + const +{ + return mCells; +} + + +Matrix4 +operator*(const Double rhs, + const Matrix4& lhs) +{ + /* Scalar multiplication is commutative. */ + return lhs * rhs; +} + #pragma mark - Rays /* diff --git a/src/basics.h b/src/basics.h index d84c197..d915b60 100644 --- a/src/basics.h +++ b/src/basics.h @@ -119,6 +119,81 @@ Vector3 LinearCombination(const Double k1, const Vector3& v1, const Double k3, const Vector3& v3); +struct Matrix4 +{ + /** Create a 4x4 zero matrix. That is, all cells are 0. */ + static Matrix4 Zero(); + + /** Create a 4x4 identity matrix. */ + static Matrix4 Identity(); + + Matrix4(); + Matrix4(const Double cells[16]); + Matrix4(const Matrix4& rhs); + + /** + * Create a 4x4 translation matrix. A translation matrix looks like this: + * + * [1 0 0 x] + * [0 1 0 y] + * [0 0 1 z] + * [0 0 0 1] + * + * @param [in] x X translation + * @param [in] y Y translation + * @param [in] z Z translation + * @returns The translation matrix + */ + static Matrix4 Translation(Double x, Double y, Double z); + + /** + * Create a 4x4 rotation matrix. A rotation matrices are quite complicated. + * + * @param [in] x X rotation angle in radians + * @param [in] y Y rotation angle in radians + * @param [in] z Z rotation angle in radians + * @returns The rotation matrix + */ + static Matrix4 Rotation(Double x, Double y, Double z); + + /** + * Get the value of the cell at (row, col). + * + * @param [in] row The row, must be less than the matrix's width + * @param [in] col The column, must be less than the matrix's height + * @returns The value of the cell at (row, col) + */ + Double& operator()(const unsigned int row, const unsigned int col); + + /** + * Scalar multiplication. + * + * @param [in] rhs The scalar factor + * @returns A copy of this matrix, multiplied by the scalar + */ + Matrix4 operator*(const Double rhs) const; + + /** + * Scalar multiplication. Multiplies this matrix by the given factor. + * + * @param [in] rhs The scalar factor + * @returns *this + */ + Matrix4& operator*=(const Double rhs); + + Matrix4 operator*(const Matrix4& rhs) const; + Matrix4& operator*=(const Matrix4& rhs); + + const Double* CArray() const; + +private: + Double mCells[16]; +}; + + +Matrix4 operator*(const Double lhs, const Matrix4& rhs); + + struct Ray { Ray(); diff --git a/test/test_basics.cc b/test/test_basics.cc index a6395b1..374a668 100644 --- a/test/test_basics.cc +++ b/test/test_basics.cc @@ -162,3 +162,112 @@ TEST_F(Vector3Test, DotProduct) { EXPECT_EQ(131.0, v1.dot(v2)); } + +#pragma mark Matrix4 Tests + +class Matrix4Test + : public ::testing::Test +{ +public: + virtual void SetUp(); + +protected: + Matrix4 m1, m2; +}; + + +void +Matrix4Test::SetUp() +{ + const Double m1Cells[] = { 1, 2, 3, 4, + 5, 6, 7, 8, + 9, 10, 11, 12, + 13, 14, 15, 16}; + m1 = Matrix4(m1Cells); + + const Double m2Cells[] = { 1, 1, 2, 3, + 5, 8, 13, 21, + 34, 55, 89, 144, + 233, 377, 610, 987}; + m2 = Matrix4(m2Cells); +} + + +TEST(Matrix4StaticTest, Zero) +{ + Matrix4 zero = Matrix4::Zero(); + for (int i = 0; i < 16; i++) { + EXPECT_EQ(zero.CArray()[i], 0.0); + } +} + + +TEST(Matrix4StaticTest, Identity) +{ + Matrix4 id = Matrix4::Identity(); + for (int i = 0; i < 4; i++) { + for (int j = 0; j < 4; j++) { + EXPECT_EQ(id.CArray()[i * 4 + j], ((i == j) ? 1.0 : 0.0)); + } + } +} + + +TEST_F(Matrix4Test, OperatorCall) +{ + for (int i = 0; i < 4; i++) { + for (int j = 0; j < 4; j++) { + EXPECT_EQ(m1(i, j), 4 * i + j + 1); + } + } +} + + +TEST_F(Matrix4Test, ScalarMultiplication) +{ + Matrix4 p1 = m1 * 2.0; + for (int i = 0; i < 4; i++) { + for (int j = 0; j < 4; j++) { + EXPECT_EQ(p1(i, j), m1(i, j) * 2.0); + } + } + + Matrix4 p2 = 2.0 * m1; + for (int i = 0; i < 4; i++) { + for (int j = 0; j < 4; j++) { + EXPECT_EQ(p2(i, j), m1(i, j) * 2.0); + } + } +} + + +TEST_F(Matrix4Test, MatrixMultiplication) +{ + const Double p1Expect[] = {1045, 1690, 2735, 4425, + 2137, 3454, 5591, 9045, + 3229, 5218, 8447, 13665, + 4321, 6982, 11303, 18285}; + Matrix4 p1 = m1 * m2; + for (int i = 0; i < 16; i++) { + EXPECT_EQ(p1.CArray()[i], p1Expect[i]); + } + + const Double p2Expect[] = { 63, 70, 77, 84, + 435, 482, 529, 576, + 2982, 3304, 3626, 3948, + 20439, 22646, 24853, 27060}; + Matrix4 p2 = m2 * m1; + for (int i = 0; i < 16; i++) { + EXPECT_EQ(p2.CArray()[i], p2Expect[i]); + } + + /* Multiplication with the identity matrix produces the same matrix. */ + Matrix4 p3 = m1 * Matrix4::Identity(); + for (int i = 0; i < 16; i++) { + EXPECT_EQ(p3.CArray()[i], m1.CArray()[i]); + } + Matrix4 p4 = Matrix4::Identity() * m1; + for (int i = 0; i < 16; i++) { + EXPECT_EQ(p4.CArray()[i], m1.CArray()[i]); + } +}