From 8d1b806ad5b2dc620aca74a88565c9cf94fdea41 Mon Sep 17 00:00:00 2001 From: Seth Parker Date: Mon, 21 Mar 2022 00:22:17 +0000 Subject: [PATCH] Port classes and features from other projects --- .clang-tidy | 2 +- .gitlab-ci.yml | 8 +- CMakeLists.txt | 9 + README.md | 57 ++++- docs/index.md | 37 +++ docs/modules.dox | 5 + docs/pages/.empty_dir | 1 + examples/CMakeLists.txt | 8 +- examples/ImageExample.cpp | 37 +++ examples/MatExample.cpp | 55 ++++ include/educelab/core.hpp | 6 + include/educelab/core/io/ImageIO.hpp | 16 ++ include/educelab/core/types/Color.hpp | 3 +- include/educelab/core/types/Image.hpp | 132 ++++++++++ include/educelab/core/types/Mat.hpp | 177 +++++++++++++ include/educelab/core/types/Signals.hpp | 4 +- include/educelab/core/types/Vec.hpp | 59 ++++- include/educelab/core/utils/Filesystem.hpp | 50 ++++ include/educelab/core/utils/LinearAlgebra.hpp | 34 +++ include/educelab/core/utils/Math.hpp | 96 +++++++ src/Image.cpp | 237 ++++++++++++++++++ src/ImageIO.cpp | 52 ++++ tests/CMakeLists.txt | 4 + tests/src/TestFilesystem.cpp | 42 ++++ tests/src/TestImage.cpp | 158 ++++++++++++ tests/src/TestLinearAlgebra.cpp | 24 ++ tests/src/TestMat.cpp | 129 ++++++++++ tests/src/TestMath.cpp | 68 +++++ tests/src/TestVec.cpp | 9 + 29 files changed, 1504 insertions(+), 15 deletions(-) create mode 100644 docs/pages/.empty_dir create mode 100644 examples/ImageExample.cpp create mode 100644 examples/MatExample.cpp create mode 100644 include/educelab/core/io/ImageIO.hpp create mode 100644 include/educelab/core/types/Image.hpp create mode 100644 include/educelab/core/types/Mat.hpp create mode 100644 include/educelab/core/utils/Filesystem.hpp create mode 100644 include/educelab/core/utils/LinearAlgebra.hpp create mode 100644 src/Image.cpp create mode 100644 src/ImageIO.cpp create mode 100644 tests/src/TestFilesystem.cpp create mode 100644 tests/src/TestImage.cpp create mode 100644 tests/src/TestLinearAlgebra.cpp create mode 100644 tests/src/TestMat.cpp diff --git a/.clang-tidy b/.clang-tidy index e2462cf..11f71e4 100644 --- a/.clang-tidy +++ b/.clang-tidy @@ -41,7 +41,7 @@ CheckOptions: - key: readability-identifier-naming.ConstexprMethodCase value: camelBack - key: readability-identifier-naming.ConstexprVariableCase - value: UPPER_CASE + value: lower_case - key: readability-identifier-naming.EnumCase value: CamelCase - key: readability-identifier-naming.EnumConstantCase diff --git a/.gitlab-ci.yml b/.gitlab-ci.yml index 0613924..1c63694 100644 --- a/.gitlab-ci.yml +++ b/.gitlab-ci.yml @@ -24,21 +24,21 @@ - *test_script ### Tests ### -test:debian:10: +test:debian:11: extends: .build_and_test stage: test needs: [] - image: volcart/vcbuilder-debian:10_v1.static + image: volcart/vcbuilder-debian:11_v1.static variables: EXTRA_CMAKE_FLAGS: "-DEDUCE_CORE_BUILD_TESTS=ON" tags: - docker -examples:debian:10: +examples:debian:11: extends: .build stage: test needs: [] - image: volcart/vcbuilder-debian:10_v1.static + image: volcart/vcbuilder-debian:11_v1.static variables: EXTRA_CMAKE_FLAGS: "-DEDUCE_CORE_BUILD_EXAMPLES=ON" tags: diff --git a/CMakeLists.txt b/CMakeLists.txt index 9fb71b3..bf92990 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -28,18 +28,27 @@ include(CheckToNumericFP) set(public_hdrs include/educelab/core.hpp include/educelab/core/Version.hpp + include/educelab/core/io/ImageIO.hpp include/educelab/core/types/Color.hpp + include/educelab/core/types/Image.hpp + include/educelab/core/types/Mat.hpp include/educelab/core/types/Mesh.hpp + include/educelab/core/types/Signals.hpp include/educelab/core/types/Uuid.hpp include/educelab/core/types/Vec.hpp + include/educelab/core/utils/Filesystem.hpp include/educelab/core/utils/Iteration.hpp + include/educelab/core/utils/LinearAlgebra.hpp include/educelab/core/utils/Math.hpp include/educelab/core/utils/String.hpp ) configure_file(src/Version.cpp.in Version.cpp) set(srcs + ${public_hdrs} ${CMAKE_CURRENT_BINARY_DIR}/Version.cpp + src/Image.cpp + src/ImageIO.cpp src/Uuid.cpp ) diff --git a/README.md b/README.md index feb0dea..3d3c63f 100644 --- a/README.md +++ b/README.md @@ -21,7 +21,7 @@ cmake -S . -B build/ -DCMAKE_BUILD_TYPE=Release ``` ## Installation -Follow the build instrcutions above, then run the following command from the root of the source directory: +Follow the build instructions above, then run the following command from the root of the source directory: ```shell # Install the library to the system @@ -36,9 +36,21 @@ following files can be installed in this way: - `utils/Iteration.hpp` - `utils/Math.hpp` - `utils/String.hpp` +- `utils/Filesystem.hpp` + - Requires: + - `utils/String.hpp` +- `utils/LinearAlgebra.hpp` + - Requires: + - `utils/Math.hpp` + - `MatrixType` and `VectorType` which implement the `Mat` and `Vec` + interfaces. +- `types/Signals.hpp` - `types/Vec.hpp` - Requires: - `utils/Math.hpp` +- `types/Mat.hpp` + - Requires: + - `types/Vec.hpp` - `types/Color.hpp` - Requires: - `types/Vec.hpp` @@ -69,6 +81,49 @@ std::cout << v0.cross(v1) << "\n"; // "[0, 0, 1]" See [examples/VecExample.cpp](examples/VecExample.cpp) for more usage examples. +### Dense 2D matrix class +```c++ +// Input point +Vec p{0, 0, 0, 1}; +std::cout << p << "\n"; // [0, 0, 0, 1] + +// Construct a translation matrix +auto M = Matrix::Eye(); +M(0, 3) = 1.f; +M(1, 3) = 2.f; +M(2, 3) = 3.f; +std::cout << "\n" << M << "\n"; // [[1, 0, 0, 1] + // [0, 1, 0, 2] + // [0, 0, 1, 3] + // [0, 0, 0, 1]] + +// Apply transform +p = translate * p; +std::cout << p << "\n"; // [1, 2, 3, 1] +``` + +See [examples/MatExample.cpp](examples/MatExample.cpp) for more usage +examples. + +### Image class +```c++ +#include + +// Construct an image +Image image(600, 800, 3, Depth::F32); + +// Fill image with a color gradient +for (const auto [y, x] : range2D(image.height(), image.width())) { + auto r = float(x) / float(image.width() - 1); + auto g = float(y) / float(image.height() - 1); + auto b = 0.25F; + image.at(y, x) = {r, g, b}; +} +``` + +See [examples/ImageExample.cpp](examples/ImageExample.cpp) for more usage +examples. + ### Iteration utilities ```c++ #include diff --git a/docs/index.md b/docs/index.md index 946ceac..a18a598 100644 --- a/docs/index.md +++ b/docs/index.md @@ -20,6 +20,43 @@ std::cout << v0 + v1 << "\n"; // "[1, 1, 0]" std::cout << v0.cross(v1) << "\n"; // "[0, 0, 1]" ``` +### Dense 2D matrix class +```{.cpp} +// Input point +Vec p{0, 0, 0, 1}; +std::cout << p << "\n"; // [0, 0, 0, 1] + +// Construct a translation matrix +auto M = Matrix::Eye(); +M(0, 3) = 1.f; +M(1, 3) = 2.f; +M(2, 3) = 3.f; +std::cout << "\n" << M << "\n"; // [[1, 0, 0, 1] + // [0, 1, 0, 2] + // [0, 0, 1, 3] + // [0, 0, 0, 1]] + +// Apply transform +p = translate * p; +std::cout << p << "\n"; // [1, 2, 3, 1] +``` + +### Image class +```{.cpp} +#include + +// Construct an image +Image image(600, 800, 3, Depth::F32); + +// Fill image with a color gradient +for (const auto [y, x] : range2D(image.height(), image.width())) { + auto r = float(x) / float(image.width() - 1); + auto g = float(y) / float(image.height() - 1); + auto b = 0.25F; + image.at(y, x) = {r, g, b}; +} +``` + ### Iteration utilities ```{.cpp} #include diff --git a/docs/modules.dox b/docs/modules.dox index d7f85e3..7417e9b 100644 --- a/docs/modules.dox +++ b/docs/modules.dox @@ -8,6 +8,11 @@ * @brief EduceLab Core library namespace */ +/** + * @namespace educelab::linalg + * @brief Linear algebra library + */ + /** * @namespace educelab::traits * @brief Class traits diff --git a/docs/pages/.empty_dir b/docs/pages/.empty_dir new file mode 100644 index 0000000..84b23ca --- /dev/null +++ b/docs/pages/.empty_dir @@ -0,0 +1 @@ +This is a placeholder file for an empty directory. diff --git a/examples/CMakeLists.txt b/examples/CMakeLists.txt index 7fcd749..af05b10 100644 --- a/examples/CMakeLists.txt +++ b/examples/CMakeLists.txt @@ -11,4 +11,10 @@ add_executable(educelab_core_UuidExample UuidExample.cpp) target_link_libraries(educelab_core_UuidExample educelab::core) add_executable(educelab_core_VecExample VecExample.cpp) -target_link_libraries(educelab_core_VecExample educelab::core) \ No newline at end of file +target_link_libraries(educelab_core_VecExample educelab::core) + +add_executable(educelab_core_ImageExample ImageExample.cpp) +target_link_libraries(educelab_core_ImageExample educelab::core) + +add_executable(educelab_core_MatExample MatExample.cpp) +target_link_libraries(educelab_core_MatExample educelab::core) \ No newline at end of file diff --git a/examples/ImageExample.cpp b/examples/ImageExample.cpp new file mode 100644 index 0000000..b6296b8 --- /dev/null +++ b/examples/ImageExample.cpp @@ -0,0 +1,37 @@ +#include "educelab/core/io/ImageIO.hpp" +#include "educelab/core/types/Image.hpp" +#include "educelab/core/types/Vec.hpp" +#include "educelab/core/utils/Iteration.hpp" + +using namespace educelab; + +auto main() -> int +{ + // Construct an image + Image image(600, 800, 3, Depth::F32); + + // Fill image with a color gradient + for (const auto [y, x] : range2D(image.height(), image.width())) { + // Create RGB values + auto r = float(x) / float(image.width() - 1); + auto g = float(y) / float(image.height() - 1); + auto b = 0.25F; + + // Assign value to the correct pixel. Like its equivalent in OpenCV, + // Image.at gives you access to the raw image buffer memory, + // reinterpreted as the type you request. Make sure you know the data + // type and channels when you access a pixel! + image.at(y, x) = {r, g, b}; + } + + // (Optional) Apply gamma correction + image = Image::Gamma(image); + + // (Optional) Convert image to the correct format for output. write_image() + // will automatically convert if the output format doesn't support the + // given bit depth, but it's always good to control the conversion. + image = image.convert(Depth::U8); + + // Write the image + write_image("educelab_core_ImageExample.ppm", image); +} diff --git a/examples/MatExample.cpp b/examples/MatExample.cpp new file mode 100644 index 0000000..c02ad94 --- /dev/null +++ b/examples/MatExample.cpp @@ -0,0 +1,55 @@ +#include +#include + +#include "educelab/core/types/Mat.hpp" +#include "educelab/core/types/Vec.hpp" +#include "educelab/core/utils/Iteration.hpp" +#include "educelab/core/utils/Math.hpp" + +using namespace educelab; + +auto main() -> int +{ + // Type aliases + using Matrix = Mat<4, 4, float>; + using Point = Vec; + + // Input point + Point p{0, 0, 0, 1}; + std::cout << "Starting point: " << p << "\n\n"; + + // Translate 1 unit along x + auto translate = Matrix::Eye(); + translate(0, 3) = 1.f; + p = translate * p; + std::cout << "Translation matrix:\n" << translate << "\n"; + std::cout << "After translation: " << p << "\n\n"; + + // Rotate around z + auto theta = to_radians(90.F); + auto rotate = Matrix::Eye(); + rotate(0, 0) = std::cos(theta); + rotate(0, 1) = -std::sin(theta); + rotate(1, 0) = std::sin(theta); + rotate(1, 1) = std::cos(theta); + p = rotate * p; + std::cout << "Rotation matrix:\n" << rotate << "\n"; + std::cout << "After rotation: " << p << "\n\n"; + + // Scale by 10 + auto scale = Matrix::Eye(); + for (auto i : range(3)) { + scale(i, i) = 10.F; + } + p = scale * p; + std::cout << "Scale matrix:\n" << scale << "\n"; + std::cout << "After scale: " << p << "\n\n"; + + // Restore original point and apply a single multiplication + p = {0, 0, 0, 1}; + auto transform = scale * rotate * translate; + std::cout << "Restored starting point: " << p << "\n"; + p = transform * p; + std::cout << "Transform matrix:\n" << transform << "\n"; + std::cout << "After transform: " << p << "\n"; +} \ No newline at end of file diff --git a/include/educelab/core.hpp b/include/educelab/core.hpp index 41f0048..3bd1c55 100644 --- a/include/educelab/core.hpp +++ b/include/educelab/core.hpp @@ -2,12 +2,18 @@ #include "educelab/core/Version.hpp" +#include "educelab/core/io/ImageIO.hpp" + #include "educelab/core/types/Color.hpp" +#include "educelab/core/types/Image.hpp" +#include "educelab/core/types/Mat.hpp" #include "educelab/core/types/Mesh.hpp" #include "educelab/core/types/Signals.hpp" #include "educelab/core/types/Uuid.hpp" #include "educelab/core/types/Vec.hpp" +#include "educelab/core/utils/Filesystem.hpp" #include "educelab/core/utils/Iteration.hpp" +#include "educelab/core/utils/LinearAlgebra.hpp" #include "educelab/core/utils/Math.hpp" #include "educelab/core/utils/String.hpp" \ No newline at end of file diff --git a/include/educelab/core/io/ImageIO.hpp b/include/educelab/core/io/ImageIO.hpp new file mode 100644 index 0000000..af6dac4 --- /dev/null +++ b/include/educelab/core/io/ImageIO.hpp @@ -0,0 +1,16 @@ +#pragma once + +/** @file */ + +#include +#include + +#include "educelab/core/types/Image.hpp" + +namespace educelab +{ + +/** @brief Write an Image to disk */ +void write_image(const std::filesystem::path& path, const Image& image); + +} // namespace educelab \ No newline at end of file diff --git a/include/educelab/core/types/Color.hpp b/include/educelab/core/types/Color.hpp index 1597359..c14a633 100644 --- a/include/educelab/core/types/Color.hpp +++ b/include/educelab/core/types/Color.hpp @@ -60,7 +60,8 @@ class Color F32C3, /** 32-bit float RGBA color */ F32C4, - /** Hexadecimal RGB color string of 3 or 6 digits (`#0a3` or `#00aa33`) + /** + * Hexadecimal RGB color string of 3 or 6 digits (`#0a3` or `#00aa33`) */ HexCode }; diff --git a/include/educelab/core/types/Image.hpp b/include/educelab/core/types/Image.hpp new file mode 100644 index 0000000..2fa4fee --- /dev/null +++ b/include/educelab/core/types/Image.hpp @@ -0,0 +1,132 @@ +#pragma once + +/** @file */ + +#include +#include + +namespace educelab +{ + +/** @brief Image bit depth */ +enum class Depth { + None, /** Unset or unspecified */ + U8, /** Unsigned 8-bit integer */ + U16, /** Unsigned 16-bit integer */ + F32 /** 32-bit float */ +}; + +/** @brief Class for storing image data */ +class Image +{ +public: + /** + * @brief Default constructor + * + * Initializes an empty image. + */ + Image() = default; + + /** + * @brief Construct a new image of given height, width, and channels + * + * All pixels are initialized to 0 in all channels. + */ + Image(std::size_t height, std::size_t width, std::size_t cns, Depth type); + + /** @brief Image width/columns */ + [[nodiscard]] auto width() const -> std::size_t; + + /** @brief Image height/rows */ + [[nodiscard]] auto height() const -> std::size_t; + + /** @brief Image aspect ratio */ + [[nodiscard]] auto aspect() const -> float; + + /** @brief Number of channels in the image */ + [[nodiscard]] auto channels() const -> std::size_t; + + /** @brief Fundamental type of each pixel element */ + [[nodiscard]] auto type() const -> Depth; + + /** @brief Whether image is empty */ + [[nodiscard]] auto empty() const -> bool; + + /** @copydoc empty() */ + explicit operator bool() const noexcept; + + /** @brief Reset image */ + void clear(); + + /** @brief Access pixel at coordinate (x,y) */ + template + auto at(std::size_t y, std::size_t x) -> T&; + + /** @copydoc at(std::size_t y, std::size_t x) */ + template + auto at(std::size_t y, std::size_t x) const -> const T&; + + /** + * @brief Convert an Image to one of the specified bit depth + * + * When converting between integer types or to a floating point image, + * intensities are scaled using the min/max ranges of the fundamental types, + * with \f$[0, 1]\f$ being the output range for a float image. When + * converting from a floating point image, the \f$[0, 1]\f$ range is scaled + * to the range of the fundamental output type. Values which would fall + * outside of the output range after scaling are clamped. + */ + static auto Convert(const Image& i, Depth type) -> Image; + + /** + * @brief Apply gamma correction to an image + * + * Applies gamma correction to each pixel in the image: + * \f$ v_{out} = v_{in}^{1/\gamma} \f$ + */ + static auto Gamma(const Image& i, float gamma = 2.F) -> Image; + + /** + * @copybrief Convert(const Image&, Depth) + * + * This is a convenience member function for calling Convert(). + */ + [[nodiscard]] auto convert(Depth type) const -> Image; + + /** @brief Returns a pointer to the underlying storage */ + [[nodiscard]] auto data() noexcept -> std::byte*; + + /** @copydoc data() */ + [[nodiscard]] auto data() const noexcept -> const std::byte*; + +private: + /** %Image height */ + std::size_t h_{0}; + /** %Image width */ + std::size_t w_{0}; + /** Number of channels */ + std::size_t cns_{0}; + /** Number of bytes between each element in a pixel */ + std::size_t stride_{0}; + /** Pixel element bit depth */ + Depth type_{Depth::None}; + /** Unravel a 2D index to a 1D position in the storage array */ + [[nodiscard]] auto unravel_(std::size_t y, std::size_t x) const + -> std::size_t; + /** Storage array */ + std::vector data_; +}; + +template +auto Image::at(std::size_t y, std::size_t x) -> T& +{ + return reinterpret_cast(data_.at(unravel_(y, x))); +} + +template +auto Image::at(std::size_t y, std::size_t x) const -> const T& +{ + return reinterpret_cast(data_.at(unravel_(y, x))); +} + +} // namespace educelab \ No newline at end of file diff --git a/include/educelab/core/types/Mat.hpp b/include/educelab/core/types/Mat.hpp new file mode 100644 index 0000000..f747035 --- /dev/null +++ b/include/educelab/core/types/Mat.hpp @@ -0,0 +1,177 @@ +#pragma once + +/** @file */ + +#include + +#include "educelab/core/types/Vec.hpp" + +namespace educelab +{ + +/** + * @brief Dense 2D matrix class for linear algebra + * + * @tparam Rows Number of rows + * @tparam Cols Number of columns + * @tparam T Fundamental element type + */ +template +class Mat +{ +public: + /** @brief Number of rows */ + static constexpr std::size_t rows{Rows}; + /** @brief Number of columns */ + static constexpr std::size_t cols{Cols}; + + /** Default constructor */ + Mat() { vals_.fill(0); }; + + /** @brief Constructor with fill values */ + template + explicit Mat(Args&&... args) + { + static_assert( + sizeof...(args) == Rows * Cols, "Incorrect number of arguments"); + std::size_t i{0}; + ((vals_[i++] = args), ...); + } + + /** @brief Matrix element access with bounds checking */ + auto at(std::size_t y, std::size_t x) -> T& + { + return vals_.at(Unravel(y, x)); + } + + /** @copydoc at(std::size_t, std::size_t) */ + auto at(std::size_t y, std::size_t x) const -> const T& + { + return vals_.at(Unravel(y, x)); + } + + /** @brief Matrix element access without bounds checking */ + auto operator()(std::size_t y, std::size_t x) -> T& + { + return vals_[Unravel(y, x)]; + } + + /** @copydoc operator()(std::size_t, std::size_t) */ + auto operator()(std::size_t y, std::size_t x) const -> const T& + { + return vals_[Unravel(y, x)]; + } + + /** @brief Return a transposed copy of the matrix */ + auto t() const -> Mat + { + Mat m; + for (std::size_t y{0}; y < Rows; y++) { + for (std::size_t x{0}; x < Cols; x++) { + m(x, y) = vals_[Unravel(y, x)]; + } + } + return m; + } + + /** @brief Construct a new identity matrix */ + static auto Eye() -> Mat + { + static_assert(Rows == Cols, "Matrix must be square"); + Mat m; + for (std::size_t i{0}; i < Rows; i++) { + m(i, i) = 1; + } + return m; + } + + /** Access to underlying storage array */ + constexpr auto data() noexcept -> T* { return vals_.data(); } + /** @copydoc data() */ + constexpr auto data() const noexcept -> const T* { return vals_.data(); } + +private: + /** Compute flat index */ + static inline auto Unravel(std::size_t y, std::size_t x) + { + return y * Cols + x; + } + /** Storage array */ + std::array vals_{}; +}; + +/** @brief Matrix-matrix multiplication operator */ +template +auto operator*(const Mat& A, const Mat& B) -> Mat +{ + Mat res; + for (std::size_t m{0}; m < M; m++) { + for (std::size_t p{0}; p < P; p++) { + for (std::size_t n{0}; n < N; n++) { + res(m, p) += A(m, n) * B(n, p); + } + } + } + return res; +} + +/** @brief Matrix-vector multiplication */ +template +auto operator*(const Mat& mat, const Vec& vec) -> Vec +{ + Vec res; + for (std::size_t m{0}; m < M; m++) { + for (std::size_t n{0}; n < N; n++) { + res[m] += mat(m, n) * vec[n]; + } + } + return res; +} + +/** @brief Calculate the 2x2 matrix determinant */ +template +auto determinant(const Mat<2, 2, T>& m) -> T +{ + return m.data()[0] * m.data()[3] - m.data()[1] * m.data()[2]; +} + +/** @brief Calculate the 3x3 matrix determinant */ +template +auto determinant(const Mat<3, 3, T>& m) -> T +{ + // a(ei − fh) − b(di − fg) + c(dh − eg) + // 0(48 - 57) - 1(38 - 56) + 2(37 - 46) + return m.data()[0] * + (m.data()[4] * m.data()[8] - m.data()[5] * m.data()[7]) - + m.data()[1] * + (m.data()[3] * m.data()[8] - m.data()[5] * m.data()[6]) + + m.data()[2] * + (m.data()[3] * m.data()[7] - m.data()[4] * m.data()[6]); +} + +/** Debug: Print a matrix to a std::ostream */ +template +auto operator<<(std::ostream& os, const Mat& mat) + -> std::ostream& +{ + os << "["; + for (std::size_t y{0}; y < Rows; y++) { + if (y != 0) { + os << " "; + } + os << "["; + for (std::size_t x{0}; x < Cols; x++) { + if (x > 0) { + os << ", "; + } + os << mat(y, x); + } + os << "]"; + if (y != Rows - 1) { + os << "\n"; + } + } + os << "]"; + return os; +} +} // namespace educelab \ No newline at end of file diff --git a/include/educelab/core/types/Signals.hpp b/include/educelab/core/types/Signals.hpp index 41bd803..3885486 100644 --- a/include/educelab/core/types/Signals.hpp +++ b/include/educelab/core/types/Signals.hpp @@ -116,7 +116,9 @@ class Signal private: /** Basic connection struct */ struct Connection { + /** Constructor */ explicit Connection(const SlotFnType& f) : slot{f} {} + /** Slot function reference */ SlotFnType slot; }; /** List of connections */ @@ -138,4 +140,4 @@ inline void Signal<>::operator()() const send(); } -} // namespace educelab \ No newline at end of file +} // namespace educelab diff --git a/include/educelab/core/types/Vec.hpp b/include/educelab/core/types/Vec.hpp index e643f50..25c31fb 100644 --- a/include/educelab/core/types/Vec.hpp +++ b/include/educelab/core/types/Vec.hpp @@ -281,7 +281,7 @@ class Vec return *this; } - /** @brief Multiplication operator */ + /** @brief Vector-scalar multiplication operator */ template < typename T2, std::enable_if_t::value, bool> = true> @@ -291,6 +291,16 @@ class Vec return lhs; } + /** @brief Scalar-vector multiplication operator */ + template < + typename T2, + std::enable_if_t::value, bool> = true> + friend auto operator*(const T2& lhs, Vec rhs) -> Vec + { + rhs *= lhs; + return rhs; + } + /** @brief Negation operator */ friend auto operator-(Vec rhs) -> Vec { @@ -310,7 +320,14 @@ class Vec return *this; } - /** @brief Division operator */ + /** + * @brief Vector-scalar division operator + * + * Performs the following operation: + * \f[ + * ret_i = \frac{rhs_i}{lhs},\ \forall i \in rhs + * \f] + */ template < typename T2, std::enable_if_t::value, bool> = true> @@ -320,9 +337,28 @@ class Vec return lhs; } + /** + * @brief Scalar-vector division operator + * + * Performs the following operation: + * \f[ + * ret_i = \frac{lhs}{rhs_i},\ \forall i \in rhs + * \f] + */ + template < + typename T2, + std::enable_if_t::value, bool> = true> + friend auto operator/(const T2& lhs, Vec rhs) -> Vec + { + for (auto& v : rhs) { + v = lhs / v; + } + return rhs; + } + /** @brief Compute the vector dot product (i.e. inner product) */ template - auto dot(const Vector& v) -> T + auto dot(const Vector& v) const -> T { return educelab::dot(val_, v); } @@ -334,14 +370,14 @@ class Vec template < typename T2, std::enable_if_t::value, bool> = true> - auto dot(const std::initializer_list& b) -> T + auto dot(const std::initializer_list& b) const -> T { return educelab::dot(val_, b); } /** @brief Compute the vector cross product */ template - auto cross(const Vector& v) -> std::enable_if_t + auto cross(const Vector& v) const -> std::enable_if_t { return educelab::cross(*this, v); } @@ -351,7 +387,7 @@ class Vec typename T2, std::size_t D = Dims, std::enable_if_t::value, bool> = true> - auto cross(const std::initializer_list& b) + auto cross(const std::initializer_list& b) const -> std::enable_if_t { return educelab::cross(*this, b); @@ -360,6 +396,15 @@ class Vec /** @brief Compute the vector magnitude */ auto magnitude() const -> T { return educelab::norm(*this, Norm::L2); } + /** @brief Compute the squared vector magnitude */ + auto magnitude2() const -> T + { + auto sum = std::accumulate( + std::begin(val_), std::end(val_), T(0), + [](auto a, auto b) { return a + (b * b); }); + return sum; + } + /** @brief Return the unit vector of this vector */ auto unit() const -> Vec { return educelab::normalize(*this); } @@ -368,6 +413,8 @@ class Vec Container val_{}; }; +/** @brief 3D, 8-bit unsigned int vector */ +using Vec3b = Vec; /** @brief 3D, 32-bit float vector */ using Vec3f = Vec; /** @brief 3D, 64-bit float vector */ diff --git a/include/educelab/core/utils/Filesystem.hpp b/include/educelab/core/utils/Filesystem.hpp new file mode 100644 index 0000000..1a10c39 --- /dev/null +++ b/include/educelab/core/utils/Filesystem.hpp @@ -0,0 +1,50 @@ +#pragma once + +/** @file */ + +#include +#include +#include +#include + +#include "educelab/core/utils/String.hpp" + +namespace educelab +{ + +namespace detail +{ +/** Compare file extension against list of extensions */ +template +auto extension_is_in_list( + std::string_view ext, const Extension& first, const ExtensionList&... rest) + -> bool +{ + if (ext == to_lower(first)) { + return true; + } + + if constexpr (sizeof...(rest) > 0) { + return extension_is_in_list(ext, rest...); + } + + return false; +} +} // namespace detail + +/** + * @brief Returns true if `path` has a file extension in `exts`. Comparison is + * case insensitive. + */ +template +auto is_file_type( + const std::filesystem::path& path, const ExtensionList&... exts) -> bool +{ + if (not path.has_extension()) { + return false; + } + auto ext = to_lower(path.extension().string().substr(1)); + return detail::extension_is_in_list(ext, exts...); +} + +} // namespace educelab diff --git a/include/educelab/core/utils/LinearAlgebra.hpp b/include/educelab/core/utils/LinearAlgebra.hpp new file mode 100644 index 0000000..274cdc1 --- /dev/null +++ b/include/educelab/core/utils/LinearAlgebra.hpp @@ -0,0 +1,34 @@ +#pragma once + +/** @file */ + +#include "educelab/core/utils/Math.hpp" + +namespace educelab::linalg +{ + +/** + * @brief Solve a system of linear equations using Cramer's rule + * + * @throws std::runtime_error if \f$ det(A) \approx 0 \f$ + */ +template +auto solve_cramer(const MatrixType& A, const VectorType& b) -> VectorType +{ + auto detA = determinant(A); + if (almost_zero(detA)) { + throw std::runtime_error("Determinant of A is zero"); + } + + VectorType res; + for (std::size_t x{0}; x < A.cols; x++) { + auto mC = A; + mC(0, x) = b[0]; + mC(1, x) = b[1]; + mC(2, x) = b[2]; + res[x] = determinant(mC) / detA; + } + return res; +} + +} // namespace educelab::linalg \ No newline at end of file diff --git a/include/educelab/core/utils/Math.hpp b/include/educelab/core/utils/Math.hpp index f193b39..293c720 100644 --- a/include/educelab/core/utils/Math.hpp +++ b/include/educelab/core/utils/Math.hpp @@ -6,6 +6,7 @@ #include #include #include +#include namespace educelab { @@ -49,6 +50,17 @@ auto cross(const T1& a, const T2& b) -> T1 return c; } +/** @brief Element-wise vector product */ +template +auto schur_product(const T1& a, const T2& b) -> T1 +{ + T1 res; + std::transform( + std::begin(a), std::end(a), std::begin(b), std::begin(res), + std::multiplies()); + return res; +} + /** @brief Vector cross product for initializer list */ template < typename T1, @@ -134,4 +146,88 @@ constexpr auto to_degrees(T2 rad) -> T return rad * T(180) / PI; } +/** + * @brief Generate a uniformly random number in the range [min, max) + * + * @note For each type T, this function statically defines a single std::mt19937 + * generator and seeds it from std::random_device. This behavior may not be + * desirable for all applications. + */ +template +inline auto random(T min = 0, T max = 1) -> T +{ + std::uniform_real_distribution distribution(min, max); + static std::random_device source; + static std::mt19937 generator(source()); + return distribution(generator); +} + +/** @brief Check if the given value is almost zero using an absolute epsilon */ +template < + typename T, + std::enable_if_t::value, bool> = true> +inline auto almost_zero(T v, T eps = 1e-7) -> bool +{ + return std::abs(v) < eps; +} + +/** + * @brief Solve for the solutions of a quadratic equation + * + * @note t0 and t1 are sorted from lowest to highest value. + */ +template < + typename T, + std::enable_if_t::value, bool> = true> +inline auto solve_quadratic(T a, T b, T c) +{ + // Zero a means the equation is linear + if (almost_zero(a)) { + throw std::invalid_argument("First quadratic coefficient is zero"); + } + + /** Structure for storing solutions to a quadratic equation */ + struct quadratic_result { + /** Whether the solution is real or imaginary */ + bool is_real{false}; + /** @copydoc is_real */ + explicit operator bool() const noexcept { return is_real; } + /** First solution */ + T t0{INF}; + /** Second solution */ + T t1{INF}; + }; + + // Calculate the discriminant + auto dis = std::pow(b, 2) - 4 * a * c; + if (dis < 0) { + return quadratic_result{}; + } + + // Setup result + quadratic_result res; + res.is_real = true; + + // Zero discriminant means t0 and t1 are the same + if (almost_zero(dis)) { + res.t0 = res.t1 = -b / (T(2) * a); + } + + // Numerically stable solution calculations + // https://pbr-book.org/3ed-2018/Utilities/Mathematical_Routines#Quadratic + else { + auto disR = std::sqrt(dis); + T q = (b < T(0)) ? T(-0.5) * (b - disR) : T(-0.5) * (b + disR); + res.t0 = q / a; + res.t1 = c / q; + } + + // Sort solutions least to greatest + if (res.t0 > res.t1) { + std::swap(res.t0, res.t1); + } + + return res; +} + } // namespace educelab \ No newline at end of file diff --git a/src/Image.cpp b/src/Image.cpp new file mode 100644 index 0000000..984fcea --- /dev/null +++ b/src/Image.cpp @@ -0,0 +1,237 @@ +#include "educelab/core/types/Image.hpp" + +#include +#include +#include +#include +#include + +using namespace educelab; + +// Conversion constants +static constexpr float max_u8{std::numeric_limits::max()}; +static constexpr float max_u16{std::numeric_limits::max()}; +static constexpr float u8_to_u16{max_u16 / max_u8}; +static constexpr float u8_to_f32{1.F / max_u8}; +static constexpr float u16_to_u8{max_u8 / max_u16}; +static constexpr float u16_to_f32{1.F / max_u16}; + +static inline auto size_of(Depth d) -> std::size_t +{ + switch (d) { + case Depth::None: + return 0; + case Depth::U8: + return 1; + case Depth::U16: + return 2; + case Depth::F32: + return 4; + } + return 0; +} + +template +void depth_cast(const std::byte* in, std::byte* out); + +template <> +void depth_cast( + const std::byte* in, std::byte* out) +{ + std::uint8_t tmpI{0}; + std::memcpy(&tmpI, in, sizeof(std::uint8_t)); + auto f = static_cast(tmpI) * u8_to_u16; + auto tmpO = static_cast(std::max(std::min(f, max_u16), 0.F)); + std::memcpy(out, &tmpO, sizeof(std::uint16_t)); +} + +template <> +void depth_cast(const std::byte* in, std::byte* out) +{ + std::uint8_t tmpI{0}; + std::memcpy(&tmpI, in, sizeof(std::uint8_t)); + auto tmpO = static_cast(tmpI) * u8_to_f32; + std::memcpy(out, &tmpO, sizeof(float)); +} + +template <> +void depth_cast( + const std::byte* in, std::byte* out) +{ + std::uint16_t tmpI{0}; + std::memcpy(&tmpI, in, sizeof(std::uint16_t)); + auto f = static_cast(tmpI) * u16_to_u8; + auto tmpO = static_cast(std::max(std::min(f, max_u8), 0.F)); + std::memcpy(out, &tmpO, sizeof(std::uint8_t)); +} + +template <> +void depth_cast(const std::byte* in, std::byte* out) +{ + std::uint16_t tmpI{0}; + std::memcpy(&tmpI, in, sizeof(std::uint16_t)); + auto tmpO = static_cast(tmpI) * u16_to_f32; + std::memcpy(out, &tmpO, sizeof(float)); +} + +template <> +void depth_cast(const std::byte* in, std::byte* out) +{ + float tmpI{0}; + std::memcpy(&tmpI, in, sizeof(float)); + auto tmpO = static_cast( + std::max(std::min(tmpI * max_u8, max_u8), 0.F)); + std::memcpy(out, &tmpO, sizeof(std::uint8_t)); +} + +template <> +void depth_cast(const std::byte* in, std::byte* out) +{ + float tmpI{0}; + std::memcpy(&tmpI, in, sizeof(float)); + auto tmpO = static_cast( + std::max(std::min(tmpI * max_u16, max_u16), 0.F)); + std::memcpy(out, &tmpO, sizeof(std::uint16_t)); +} + +template +void pixel_cast(const std::byte* in, std::byte* out, std::size_t cns) +{ + for (std::size_t idx{0}; idx < cns; idx++) { + auto i_idx = idx * sizeof(TIn); + auto o_idx = idx * sizeof(TOut); + depth_cast(&in[i_idx], &out[o_idx]); + } +} + +static void convert_pixel( + const std::byte* in, + Depth inType, + std::byte* out, + Depth outType, + std::size_t cns) +{ + + // 8bpc -> 16bpc + if (inType == Depth::U8 and outType == Depth::U16) { + pixel_cast(in, out, cns); + } + + // 8bpc -> float + else if (inType == Depth::U8 and outType == Depth::F32) { + pixel_cast(in, out, cns); + } + + // 16bpc -> 8bpc + else if (inType == Depth::U16 and outType == Depth::U8) { + pixel_cast(in, out, cns); + } + + // 16bpc -> float + else if (inType == Depth::U16 and outType == Depth::F32) { + pixel_cast(in, out, cns); + } + + // float -> 8bpc + else if (inType == Depth::F32 and outType == Depth::U8) { + pixel_cast(in, out, cns); + } + + // float -> 16bpc + else if (inType == Depth::F32 and outType == Depth::U16) { + pixel_cast(in, out, cns); + } + + // error + else { + throw std::runtime_error("Conversion not supported."); + } +} + +Image::Image(std::size_t height, std::size_t width, std::size_t cns, Depth type) + : h_{height} + , w_{width} + , cns_{cns} + , stride_{size_of(type)} + , type_{type} + , data_(h_ * w_ * cns_ * stride_, std::byte{0}) +{ + assert(h_ > 0 and w_ > 0 and cns_ > 0); +} + +auto Image::width() const -> std::size_t { return w_; } + +auto Image::height() const -> std::size_t { return h_; } + +auto Image::aspect() const -> float +{ + return (h_ == 0) ? 0 : static_cast(w_) / static_cast(h_); +} + +auto Image::channels() const -> std::size_t { return cns_; } + +auto Image::type() const -> Depth { return type_; } + +auto Image::empty() const -> bool { return data_.empty(); } + +Image::operator bool() const noexcept { return not empty(); } + +void Image::clear() +{ + h_ = 0; + w_ = 0; + cns_ = 0; + stride_ = 0; + type_ = Depth::None; + data_.clear(); +} + +auto Image::Convert(const Image& i, Depth type) -> Image +{ + // Return if image is already of requested type + if (i.type() == type) { + return i; + } + + // Setup output image + Image result(i.h_, i.w_, i.cns_, type); + + // Iterate over original image data + for (std::size_t y{0}; y < i.h_; y++) { + for (std::size_t x{0}; x < i.w_; x++) { + auto i_idx = i.unravel_(y, x); + auto r_idx = result.unravel_(y, x); + convert_pixel( + &i.data_[i_idx], i.type(), &result.data_[r_idx], type, i.cns_); + } + } + + return result; +} + +auto Image::Gamma(const Image& i, float gamma) -> Image +{ + auto result = Convert(i, Depth::F32); + for (std::size_t y{0}; y < result.h_; y++) { + for (std::size_t x{0}; x < result.w_; x++) { + auto idx = result.unravel_(y, x); + for (std::size_t c{0}; c < result.cns_; c++) { + auto o = c * sizeof(float); + auto* v = reinterpret_cast(&result.data_[idx + o]); + *v = std::pow(*v, 1.F / gamma); + } + } + } + return Convert(result, i.type()); +} + +auto Image::convert(Depth type) const -> Image { return Convert(*this, type); } + +auto Image::data() noexcept -> std::byte* { return data_.data(); } + +auto Image::data() const noexcept -> const std::byte* { return data_.data(); } + +auto Image::unravel_(std::size_t y, std::size_t x) const -> std::size_t +{ + return (y % h_) * w_ * cns_ * stride_ + x * cns_ * stride_; +} diff --git a/src/ImageIO.cpp b/src/ImageIO.cpp new file mode 100644 index 0000000..2c12c76 --- /dev/null +++ b/src/ImageIO.cpp @@ -0,0 +1,52 @@ +#include "educelab/core/io/ImageIO.hpp" + +#include "educelab/core/types/Vec.hpp" +#include "educelab/core/utils/Filesystem.hpp" +#include "educelab/core/utils/Iteration.hpp" + +using namespace educelab; +namespace el = educelab; +namespace fs = std::filesystem; + +// Write PPM image +static void ppm_write(const fs::path& path, const Image& image) +{ + // PPM only supports 8bpc, RGB + if (image.channels() != 3) { + throw std::runtime_error( + "Unsupported number of channels: " + + std::to_string(image.channels())); + } + auto tmp = image.convert(Depth::U8); + + // Open output file + std::ofstream file(path); + if (not file.is_open()) { + throw std::runtime_error("Cannot open file: " + path.string()); + } + + // PPM Header + file << "P3\n"; + file << tmp.width() << ' ' << tmp.height() << '\n'; + file << "255\n"; + + // Pixels + for (const auto [y, x] : range2D(tmp.height(), tmp.width())) { + auto p = tmp.at(y, x); + file << static_cast(p[0]) << ' '; + file << static_cast(p[1]) << ' '; + file << static_cast(p[2]) << '\n'; + } + + file.close(); +} + +void el::write_image(const fs::path& path, const Image& image) +{ + if (is_file_type(path, "ppm")) { + ppm_write(path, image); + } else { + auto ext = path.extension().string(); + throw std::invalid_argument("Unsupported file type: " + ext); + } +} diff --git a/tests/CMakeLists.txt b/tests/CMakeLists.txt index 55834a5..bb16868 100644 --- a/tests/CMakeLists.txt +++ b/tests/CMakeLists.txt @@ -18,7 +18,11 @@ endif() ## Build the tests ## set(tests src/TestColor.cpp + src/TestFilesystem.cpp + src/TestImage.cpp src/TestIteration.cpp + src/TestLinearAlgebra.cpp + src/TestMat.cpp src/TestMath.cpp src/TestMesh.cpp src/TestSignals.cpp diff --git a/tests/src/TestFilesystem.cpp b/tests/src/TestFilesystem.cpp new file mode 100644 index 0000000..685414b --- /dev/null +++ b/tests/src/TestFilesystem.cpp @@ -0,0 +1,42 @@ +#include + +#include "educelab/core/utils/Filesystem.hpp" + +using namespace educelab; +namespace fs = std::filesystem; + +TEST(Filesystem, IsFileTypeString) +{ + // Single matching extension + EXPECT_TRUE(is_file_type("some/path.jpg", "jpg")); + // Double matching extension + EXPECT_TRUE(is_file_type("some/path.jpg", "tif", "jpg")); + + // Single non-matching extension + EXPECT_FALSE(is_file_type("some/path.jpg", "tif")); + // Double non-matching extension + EXPECT_FALSE(is_file_type("some/path.jpg", "tif", "bmp")); + + // No extension + EXPECT_FALSE(is_file_type("file", "")); +} + +TEST(Filesystem, IsFileTypePath) +{ + // Path + fs::path input{"some/path.jpg"}; + + // Single matching extension + EXPECT_TRUE(is_file_type(input, "jpg")); + // Double matching extension + EXPECT_TRUE(is_file_type(input, "tif", "jpg")); + + // Single non-matching extension + EXPECT_FALSE(is_file_type(input, "tif")); + // Double non-matching extension + EXPECT_FALSE(is_file_type(input, "tif", "bmp")); + + // No extension + input = "file"; + EXPECT_FALSE(is_file_type(input, "")); +} \ No newline at end of file diff --git a/tests/src/TestImage.cpp b/tests/src/TestImage.cpp new file mode 100644 index 0000000..58cfb13 --- /dev/null +++ b/tests/src/TestImage.cpp @@ -0,0 +1,158 @@ +#include + +#include "educelab/core/types/Image.hpp" +#include "educelab/core/utils/Iteration.hpp" + +using namespace educelab; + +TEST(Image, DefaultConstructor) +{ + Image img; + EXPECT_EQ(img.height(), 0); + EXPECT_EQ(img.width(), 0); + EXPECT_EQ(img.channels(), 0); + EXPECT_EQ(img.type(), Depth::None); + EXPECT_FALSE(img); +} + +TEST(Image, PropertiesConstructorU8) +{ + // Construct image + Image img(5, 10, 1, Depth::U8); + EXPECT_EQ(img.height(), 5); + EXPECT_EQ(img.width(), 10); + EXPECT_EQ(img.channels(), 1); + EXPECT_EQ(img.type(), Depth::U8); + EXPECT_TRUE(img); + + // Check zero initialization + for (const auto [y, x] : range2D(img.height(), img.width())) { + EXPECT_EQ(img.at(y, x), std::uint8_t{0}); + } +} + +TEST(Image, PropertiesConstructorU16) +{ + // Construct image + Image img(5, 10, 1, Depth::U16); + EXPECT_EQ(img.height(), 5); + EXPECT_EQ(img.width(), 10); + EXPECT_EQ(img.channels(), 1); + EXPECT_EQ(img.type(), Depth::U16); + EXPECT_TRUE(img); + + // Check zero initialization + for (const auto [y, x] : range2D(img.height(), img.width())) { + EXPECT_EQ(img.at(y, x), std::uint16_t{0}); + } +} + +TEST(Image, PropertiesConstructorF32) +{ + // Construct image + Image img(5, 10, 1, Depth::F32); + EXPECT_EQ(img.height(), 5); + EXPECT_EQ(img.width(), 10); + EXPECT_EQ(img.channels(), 1); + EXPECT_EQ(img.type(), Depth::F32); + EXPECT_TRUE(img); + + // Check zero initialization + for (const auto [y, x] : range2D(img.height(), img.width())) { + EXPECT_FLOAT_EQ(img.at(y, x), 0.F); + } +} + +TEST(Image, ConvertFromU8) +{ + // Construct image + Image img(10, 10, 1, Depth::U8); + for (const auto i : range(10)) { + img.at(i, i) = 255; + img.at(i, 9 - i) = 127; + } + + // Convert to 16-bpc + auto img16 = Image::Convert(img, Depth::U16); + for (const auto i : range(10)) { + EXPECT_EQ(img16.at(i, i), 65535); + EXPECT_EQ(img16.at(i, 9 - i), 32639); + } + + // Convert to 32-bpc + auto img32 = Image::Convert(img, Depth::F32); + for (const auto i : range(10)) { + EXPECT_FLOAT_EQ(img32.at(i, i), 1.F); + EXPECT_FLOAT_EQ(img32.at(i, 9 - i), 127.F / 255.F); + } +} + +TEST(Image, ConvertFromU16) +{ + // Construct image + Image img(10, 10, 1, Depth::U16); + for (const auto i : range(10)) { + img.at(i, i) = 65535; + img.at(i, 9 - i) = 32767; + } + + // Convert to 8-bpc + auto img8 = Image::Convert(img, Depth::U8); + for (const auto i : range(10)) { + EXPECT_EQ(img8.at(i, i), 255); + EXPECT_EQ(img8.at(i, 9 - i), 127); + } + + // Convert to 32-bpc + auto img32 = Image::Convert(img, Depth::F32); + for (const auto i : range(10)) { + EXPECT_FLOAT_EQ(img32.at(i, i), 1.F); + EXPECT_FLOAT_EQ(img32.at(i, 9 - i), 32767.F / 65535.F); + } +} + +TEST(Image, ConvertFromF32) +{ + // Construct image + Image img(10, 10, 1, Depth::F32); + for (const auto i : range(10)) { + img.at(i, i) = 1.F; + img.at(i, 9 - i) = 0.5F; + } + + // Convert to 8-bpc + auto img8 = Image::Convert(img, Depth::U8); + for (const auto i : range(10)) { + EXPECT_EQ(img8.at(i, i), 255); + EXPECT_EQ(img8.at(i, 9 - i), 127); + } + + // Convert to 16-bpc + auto img16 = Image::Convert(img, Depth::U16); + for (const auto i : range(10)) { + EXPECT_EQ(img16.at(i, i), 65535); + EXPECT_EQ(img16.at(i, 9 - i), 32767); + } +} + +TEST(Image, Gamma) +{ + // Construct image and expected array + Image img(1, 11, 1, Depth::F32); + std::array expected; + + // Fill image and expected array + float gamma{2.F}; + for (const auto x : range(11)) { + // Assign to image + img.at(0, x) = 0.1F * float(x); + // Calculate expected + expected.at(x) = std::pow(0.1F * float(x), 1.F / gamma); + } + + // Apply gamma correction + auto gammaImg = Image::Gamma(img, 2.F); + for (const auto x : range(11)) { + EXPECT_FLOAT_EQ(gammaImg.at(0, x), expected.at(x)); + } +} \ No newline at end of file diff --git a/tests/src/TestLinearAlgebra.cpp b/tests/src/TestLinearAlgebra.cpp new file mode 100644 index 0000000..90874be --- /dev/null +++ b/tests/src/TestLinearAlgebra.cpp @@ -0,0 +1,24 @@ +#include + +#include "educelab/core/types/Mat.hpp" +#include "educelab/core/types/Vec.hpp" +#include "educelab/core/utils/LinearAlgebra.hpp" + +using namespace educelab; +using namespace educelab::linalg; + +TEST(LinearAlgebra, SolveCramer) +{ + Mat<3, 3> A{2, 1, 1, 1, -1, -1, 1, 2, 1}; + Vec3f b{3, 0, 0}; + Vec3f x; + EXPECT_NO_THROW(x = solve_cramer(A, b)); + EXPECT_EQ(x, Vec3f(1, -2, 3)); +} + +TEST(LinearAlgebra, SolveCramerError) +{ + Mat<3, 3> A{1, 1, 1, 1, 1, 2, 1, 1, 3}; + Vec3f b{1, 3, -1}; + EXPECT_THROW(solve_cramer(A, b), std::runtime_error); +} diff --git a/tests/src/TestMat.cpp b/tests/src/TestMat.cpp new file mode 100644 index 0000000..ce5620b --- /dev/null +++ b/tests/src/TestMat.cpp @@ -0,0 +1,129 @@ +#include + +#include "educelab/core/types/Mat.hpp" +#include "educelab/core/types/Vec.hpp" +#include "educelab/core/utils/Iteration.hpp" + +using namespace educelab; + +using Mat3f = Mat<3, 3>; + +TEST(Mat, DefaultConstructor) +{ + Mat3f m; + for (auto [y, x] : range2D(3, 3)) { + EXPECT_FLOAT_EQ(m(y, x), 0.F); + } +} + +TEST(Mat, FillConstructor) +{ + Mat3f m{0, 1, 2, 3, 4, 5, 6, 7, 8}; + for (auto i : range(9)) { + EXPECT_FLOAT_EQ(m.data()[i], float(i)); + } +} + +TEST(Mat, At) +{ + Mat3f m; + // In bounds access (non-const) + for (auto [y, x] : range2D(3, 3)) { + EXPECT_NO_THROW(m.at(y, x) = float(y * m.cols + x)); + } + + // In bounds access (const) + for (auto [y, x] : range2D(3, 3)) { + EXPECT_FLOAT_EQ(m.at(y, x), float(y * m.cols + x)); + } + + // Out of bounds access + EXPECT_THROW(m.at(3, 3), std::out_of_range); +} + +TEST(Mat, AccessOperator) +{ + Mat3f m; + // In bounds access (non-const) + for (auto [y, x] : range2D(3, 3)) { + EXPECT_NO_FATAL_FAILURE(m(y, x) = float(y * m.cols + x)); + } + + // In bounds access (const) + for (auto [y, x] : range2D(3, 3)) { + EXPECT_FLOAT_EQ(m(y, x), float(y * m.cols + x)); + } +} + +TEST(Mat, Transpose) +{ + // Build Mat + Mat3f m{0, 1, 2, 3, 4, 5, 6, 7, 8}; + m = m.t(); + + // Iterate, but transpose indices + for (auto [y, x] : range2D(3, 3)) { + EXPECT_FLOAT_EQ(m.at(x, y), float(y * m.cols + x)); + } +} + +TEST(Mat, Eye) +{ + auto m = Mat3f::Eye(); + for (auto [y, x] : range2D(3, 3)) { + EXPECT_FLOAT_EQ(m.at(y, x), (y == x) ? 1.F : 0.F); + } +} + +TEST(Mat, MatrixMatrixMultiplication) +{ + // Square multiplication + Mat<2, 2> m0{1, 2, 3, 4}; + Mat<2, 2> m1{5, 6, 7, 8}; + auto result = m0 * m1; + + // Check result + Mat<2, 2> expected{19, 22, 43, 50}; + for (auto [y, x] : range2D(2, 2)) { + EXPECT_FLOAT_EQ(result.at(y, x), expected.at(y, x)); + } + + // Non-square multiplication + Mat<2, 3> m2{1, 2, 3, 4, 5, 6}; + Mat<3, 2> m3{7, 8, 9, 10, 11, 12}; + result = m2 * m3; + + // Check result + expected = Mat<2, 2>{58, 64, 139, 154}; + for (auto [y, x] : range2D(2, 2)) { + EXPECT_FLOAT_EQ(result.at(y, x), expected.at(y, x)); + } +} + +TEST(Mat, MatrixVectorMultiplication) +{ + using Vec4f = Vec; + + // Input point + Vec4f x{0, 0, 0, 1}; + + // Translation matrix + auto M = Mat<4, 4>::Eye(); + M(0, 3) = 1.F; + M(1, 3) = 2.F; + M(2, 3) = 3.F; + auto result = M * x; + EXPECT_EQ(result, Vec4f(1, 2, 3, 1)); +} + +TEST(Mat, Determinant2x2) +{ + Mat<2, 2> m{1, 2, 3, 4}; + EXPECT_FLOAT_EQ(determinant(m), -2.F); +} + +TEST(Mat, Determinant3x3) +{ + Mat3f m{1, 2, 3, 4, 5, 6, 7, 8, 9}; + EXPECT_FLOAT_EQ(determinant(m), 0.F); +} \ No newline at end of file diff --git a/tests/src/TestMath.cpp b/tests/src/TestMath.cpp index aa00e13..201ab8b 100644 --- a/tests/src/TestMath.cpp +++ b/tests/src/TestMath.cpp @@ -65,4 +65,72 @@ TEST(Math, InteriorAngle) using Vec2f = Vec; EXPECT_EQ(interior_angle(Vec2f{1, 0}, Vec2f{0, 1}), to_radians(90)); EXPECT_EQ(interior_angle(Vec3f{1, 0, 0}, Vec3f{0, 1, 0}), to_radians(90)); +} + +TEST(Math, RandomFloatDefault) +{ + for ([[maybe_unused]] const auto r : range(1000)) { + auto val = random(); + EXPECT_GE(val, 0.F); + EXPECT_LT(val, 1.F); + } +} + +TEST(Math, RandomFloatCustom) +{ + for ([[maybe_unused]] const auto r : range(1000)) { + auto val = random(0.F, 10.F); + EXPECT_GE(val, 0.F); + EXPECT_LT(val, 10.F); + } +} + +TEST(Math, RandomDoubleDefault) +{ + for ([[maybe_unused]] const auto r : range(1000)) { + auto val = random(); + EXPECT_GE(val, 0.); + EXPECT_LT(val, 1.); + } +} + +TEST(Math, RandomDoubleCustom) +{ + for ([[maybe_unused]] const auto r : range(1000)) { + auto val = random(0., 10.); + EXPECT_GE(val, 0.); + EXPECT_LT(val, 10.); + } +} + +TEST(Math, AlmostZero) +{ + EXPECT_TRUE(almost_zero(1e-8F)); + EXPECT_FALSE(almost_zero(1e-7F)); +} + +TEST(Math, SolveQuadraticReal) +{ + if (auto res = solve_quadratic(5.F, 6.F, 1.F)) { + EXPECT_FLOAT_EQ(res.t0, -1.F); + EXPECT_FLOAT_EQ(res.t1, -0.2F); + } +} + +TEST(Math, SolveQuadraticImaginary) +{ + EXPECT_FALSE(solve_quadratic(5.F, 2.F, 1.F)); +} + +TEST(Math, SolveQuadraticLinear) +{ + EXPECT_THROW(solve_quadratic(0.F, 2.F, 1.F), std::invalid_argument); +} + +TEST(Math, SchurProduct) +{ + Vec3f a{1, 2, 3}; + Vec3f b{4, 5, 6}; + auto result = schur_product(a, b); + EXPECT_EQ(result, Vec3f(4, 10, 18)); } \ No newline at end of file diff --git a/tests/src/TestVec.cpp b/tests/src/TestVec.cpp index 591d3d4..85695a0 100644 --- a/tests/src/TestVec.cpp +++ b/tests/src/TestVec.cpp @@ -38,6 +38,7 @@ TEST(Vec, OperatorMultiply) { Vec3f a{1, 1, 1}; EXPECT_EQ(a * 2, Vec3f(2, 2, 2)); + EXPECT_EQ(2 * a, Vec3f(2, 2, 2)); EXPECT_EQ(a, Vec3f(1, 1, 1)); EXPECT_EQ(a *= 2, Vec3f(2, 2, 2)); EXPECT_EQ(a, Vec3f(2, 2, 2)); @@ -49,6 +50,7 @@ TEST(Vec, OperatorDivide) { Vec3f a{2, 2, 2}; EXPECT_EQ(a / 2, Vec3f(1, 1, 1)); + EXPECT_EQ(2 / a, Vec3f(1, 1, 1)); EXPECT_EQ(a, Vec3f(2, 2, 2)); EXPECT_EQ(a /= 2, Vec3f(1, 1, 1)); EXPECT_EQ(a, Vec3f(1, 1, 1)); @@ -89,6 +91,13 @@ TEST(Vec, Magnitude) EXPECT_EQ(Vec3f(0, 0, 3).magnitude(), 3.f); } +TEST(Vec, Magnitude2) +{ + EXPECT_EQ(Vec3f(1, 0, 0).magnitude2(), 1.f); + EXPECT_EQ(Vec3f(0, 2, 0).magnitude2(), 4.f); + EXPECT_EQ(Vec3f(0, 0, 3).magnitude2(), 9.f); +} + TEST(Vec, UnitVector) { Vec3f a{2, 0, 0};