diff --git a/common/qvec.cc b/common/qvec.cc index e5b89cb1..084298fb 100644 --- a/common/qvec.cc +++ b/common/qvec.cc @@ -137,6 +137,27 @@ qmat4x4f qv::inverse(const qmat4x4f &input) return qmat4x4f(qv::inverse(qmat4x4d(input))); } +qmat3x3f qv::inverse(const qmat3x3f &m) +{ + qmat4x4d temp {}; + for (int r = 0; r < 3; ++r) { + for (int c = 0; c < 3; ++c) { + temp.at(r,c) = m.at(r,c); + } + } + + temp = qv::inverse(temp); + + // return the upper-left 3x3 matrix + qmat3x3f result; + for (int r = 0; r < 3; ++r) { + for (int c = 0; c < 3; ++c) { + result.at(r,c) = temp.at(r,c); + } + } + return result; +} + qmat2x2f qv::inverse(const qmat2x2f &m) { // http://www.mathwords.com/i/inverse_of_a_matrix.htm diff --git a/include/common/qvec.hh b/include/common/qvec.hh index 0a16b8eb..55c5c7a1 100644 --- a/include/common/qvec.hh +++ b/include/common/qvec.hh @@ -909,6 +909,19 @@ public: std::copy(list.begin(), list.end(), m_values.begin()); } + // static factory, row-major order + static qmat row_major(std::initializer_list list) { + assert(list.size() == NRow * NCol); + + qmat result; + for (size_t i = 0; i < NRow; i++) { // for each row + for (size_t j = 0; j < NCol; j++) { // for each col + result.at(i, j) = *(list.begin() + (i*NCol + j)); + } + } + return result; + } + // Sort support [[nodiscard]] constexpr bool operator<(const qmat &other) const { return m_values < other.m_values; } [[nodiscard]] constexpr bool operator<=(const qmat &other) const { return m_values <= other.m_values; } @@ -1043,6 +1056,7 @@ namespace qv [[nodiscard]] qmat4x4d inverse(const qmat4x4d &input); [[nodiscard]] qmat2x2f inverse(const qmat2x2f &input); +[[nodiscard]] qmat3x3f inverse(const qmat3x3f &input); }; // namespace qv // returns the normalized direction from `start` to `stop` in the `dir` param diff --git a/tests/test_light.cc b/tests/test_light.cc index f256a4e3..710b79d8 100644 --- a/tests/test_light.cc +++ b/tests/test_light.cc @@ -871,6 +871,32 @@ TEST_CASE("matrix2x2inv") REQUIRE(std::isnan(nanMat.at(0, 0))); } +TEST_CASE("matrix3x3inv") +{ + std::mt19937 engine(0); + std::uniform_real_distribution dis(-4096, 4096); + + qmat3x3f randMat; + for (int i = 0; i < 3; i++) + for (int j = 0; j < 3; j++) + randMat.at(i, j) = dis(engine); + + qmat3x3f randInv = qv::inverse(randMat); + REQUIRE_FALSE(std::isnan(randInv.at(0, 0))); + + qmat3x3f prod = randMat * randInv; + for (int i = 0; i < 3; i++) { + for (int j = 0; j < 3; j++) { + float exp = (i == j) ? 1.0f : 0.0f; + REQUIRE(fabs(exp - prod.at(i, j)) < 0.001); + } + } + + // check non-invertible gives nan + qmat3x3f nanMat = qv::inverse(qmat3x3f(0)); + REQUIRE(std::isnan(nanMat.at(0, 0))); +} + TEST_CASE("matrix4x4inv") { std::mt19937 engine(0); @@ -897,6 +923,19 @@ TEST_CASE("matrix4x4inv") REQUIRE(std::isnan(nanMat.at(0, 0))); } +TEST_CASE("qmat_construct_initialize") { + const qmat2x2f test{1,2,3,4}; // column major + + CHECK(qvec2f{1,3} == test.row(0)); + CHECK(qvec2f{2,4} == test.row(1)); +} + +TEST_CASE("qmat_construct_row_major") { + const qmat2x2f test = qmat2x2f::row_major({1, 2, 3, 4}); + + CHECK(qvec2f{1,2} == test.row(0)); + CHECK(qvec2f{3,4} == test.row(1)); +} } TEST_SUITE("trace") {