diff --git a/examples/simple.cpp b/examples/simple.cpp index d7828c09..c136b44b 100644 --- a/examples/simple.cpp +++ b/examples/simple.cpp @@ -1,7 +1,6 @@ #include //std::array #include #include // std::printf -#include // std::memset #include // std::unique_ptr #define TPH_POISSON_IMPLEMENTATION diff --git a/include/thinks/tph_poisson.h b/include/thinks/tph_poisson.h index e5fd4a1d..aa36fbd3 100644 --- a/include/thinks/tph_poisson.h +++ b/include/thinks/tph_poisson.h @@ -12,7 +12,6 @@ /* * TODOS: - * - Move implementation macros to implementation section?? * - Build and run tests with sanitizers!! */ @@ -96,19 +95,44 @@ struct tph_poisson_sampling_ /* clang-format on */ /** - * @brief ... + * Generates a list of samples with the guarantees: (1) No two samples are closer to each other than + * args.radius; (2) No sample is outside the region [args.bounds_min, args.bounds_max]. + * + * The algorithm tries to fit as many samples as possible into the region without violating the + * above requirements. The samples are accessed using the tph_poisson_get_samples function. + + * If the arguments are invalid TPH_POISSON_INVALID_ARGS is returned. + * The arguments are invalid if: + * - args.radius is <= 0, or + * - args.ndims is < 1, or + * - args.bounds_min[i] >= args.bounds_max[i], or + * - args.max_sample_attempts == 0, or + * - an invalid allocator is provided. + * + * If a memory allocation fails TPH_POISSON_BAD_ALLOC is returned. + * If the number of samples exceeds the maximum number TPH_POISSON_OVERFLOW is returned. + * + * @param sampling Sampling to store samples. + * @param args Arguments. + * @param alloc Optional custom allocator (may be null). + * @return TPH_POISSON_SUCCESS if no errors; otherwise a non-zero error code. */ extern int tph_poisson_create(tph_poisson_sampling *sampling, const tph_poisson_args *args, tph_poisson_allocator *alloc); /** - * @brief ... + * @brief Frees all memory used by the sampling. + * @param sampling Sampling to store samples. */ extern void tph_poisson_destroy(tph_poisson_sampling *sampling); /** - * @brief ... + * Returns a pointer to the samples in the provided sampling. Samples are stored as N-dimensional + * points, i.e. the first N values are the coordinates of the first point, etc. Note that + * sampling.ndims and sampling.nsamples can be used to unpack the raw samples into points. + * @param sampling Sampling to store samples. + * @return Pointer to samples. */ extern const tph_poisson_real *tph_poisson_get_samples(tph_poisson_sampling *sampling); @@ -127,6 +151,12 @@ extern const tph_poisson_real *tph_poisson_get_samples(tph_poisson_sampling *sam #include // alignof +#if defined(_MSC_VER) && !defined(__cplusplus) +#define TPH_POISSON_INLINE __inline +#else +#define TPH_POISSON_INLINE inline +#endif + #ifndef tph_poisson_assert #include #define tph_poisson_assert(_X_) assert((_X_)) @@ -142,12 +172,6 @@ extern const tph_poisson_real *tph_poisson_get_samples(tph_poisson_sampling *sam #define TPH_POISSON_MEMSET(s, c, n) memset((s), (c), (n)) #endif -#if defined(_MSC_VER) && !defined(__cplusplus) -#define TPH_POISSON_INLINE __inline -#else -#define TPH_POISSON_INLINE inline -#endif - /* * MEMORY */ @@ -299,7 +323,8 @@ typedef struct tph_poisson_vec_ } tph_poisson_vec; /** - * @brief Returns non-zero if the vector is in a valid state; otherwise zero. Used for debugging. + * @brief Returns non-zero if the vector is in a valid state; otherwise zero. Only used for + * debugging. * @param ElemT Vector element type. * @param vec Vector. * @return Non-zero if the vector is in a valid state; otherwise zero. @@ -850,13 +875,14 @@ static int tph_poisson_existing_sample_within_radius(tph_poisson_context *ctx, const tph_poisson_real *cell_sample = NULL; int32_t i = -1; ptrdiff_t k = -1; + const int32_t ndims = ctx->ndims; TPH_POISSON_MEMCPY( - ctx->grid_index, min_grid_index, (size_t)(ctx->ndims * (ptrdiff_t)sizeof(ptrdiff_t))); + ctx->grid_index, min_grid_index, (size_t)(ndims * (ptrdiff_t)sizeof(ptrdiff_t))); do { /* Compute linear grid index. */ tph_poisson_assert((0 <= ctx->grid_index[0]) & (ctx->grid_index[0] < ctx->grid_size[0])); k = ctx->grid_index[0]; - for (i = 1; i < ctx->ndims; ++i) { + for (i = 1; i < ndims; ++i) { /* Not checking for overflow! */ tph_poisson_assert((0 <= ctx->grid_index[i]) & (ctx->grid_index[i] < ctx->grid_size[i])); k += ctx->grid_index[i] * ctx->grid_stride[i]; @@ -867,10 +893,10 @@ static int tph_poisson_existing_sample_within_radius(tph_poisson_context *ctx, /* Compute (squared) distance to the existing sample and then check if the existing sample is * closer than (squared) radius to the provided sample. */ cell_sample = - (const tph_poisson_real *)samples->begin + (ptrdiff_t)ctx->grid_cells[k] * ctx->ndims; + (const tph_poisson_real *)samples->begin + (ptrdiff_t)ctx->grid_cells[k] * ndims; di = sample[0] - cell_sample[0]; d_sqr = di * di; - for (i = 1; i < ctx->ndims; ++i) { + for (i = 1; i < ndims; ++i) { di = sample[i] - cell_sample[i]; d_sqr += di * di; } @@ -880,7 +906,7 @@ static int tph_poisson_existing_sample_within_radius(tph_poisson_context *ctx, /* Iterate over grid index range. Enumerate every grid index between min_grid_index and * max_grid_index (inclusive) exactly once. Assumes that min_index is element-wise less than or * equal to max_index. */ - for (i = 0; i < ctx->ndims; ++i) { + for (i = 0; i < ndims; ++i) { tph_poisson_assert(min_grid_index[i] <= max_grid_index[i]); ctx->grid_index[i]++; if (ctx->grid_index[i] <= max_grid_index[i]) { break; } @@ -889,7 +915,7 @@ static int tph_poisson_existing_sample_within_radius(tph_poisson_context *ctx, /* If the above loop ran to completion, without triggering the break, the grid_index has been * set to its original value (min_grid_index). Since this was the starting value for grid_index * we exit the outer loop when this happens. */ - } while (i != ctx->ndims); + } while (i != ndims); /* No existing sample was found to be closer to the provided sample than the radius. */ return 0; @@ -1022,17 +1048,19 @@ int tph_poisson_create(tph_poisson_sampling *s, void tph_poisson_destroy(tph_poisson_sampling *s) { - /* Sanity-check that the sampling was initialized by a call to tph_poisson_create. If the sampling - * was default-initialized do nothing. */ - if (s && s->internal) { - tph_poisson_sampling_internal *internal = s->internal; - tph_poisson_vec_free(&internal->samples, &internal->alloc); - tph_poisson_free_fn free_fn = internal->alloc.free; - void *alloc_ctx = internal->alloc.ctx; - free_fn(internal->mem, internal->mem_size, alloc_ctx); - - /* Protects from destroy being called more than once causing a double-free error. */ - s->internal = NULL; + if (s) { + s->ndims = 0; + s->nsamples = 0; + if (s->internal) { + tph_poisson_sampling_internal *internal = s->internal; + tph_poisson_vec_free(&internal->samples, &internal->alloc); + tph_poisson_free_fn free_fn = internal->alloc.free; + void *alloc_ctx = internal->alloc.ctx; + free_fn(internal->mem, internal->mem_size, alloc_ctx); + + /* Protects from destroy being called more than once causing a double-free error. */ + s->internal = NULL; + } } } @@ -1043,20 +1071,18 @@ const tph_poisson_real *tph_poisson_get_samples(tph_poisson_sampling *s) } /* Clean up internal macros. */ +#undef TPH_POISSON_INLINE #undef tph_poisson_assert #undef TPH_POISSON_MEMCPY #undef TPH_POISSON_MEMSET -#undef TPH_POISSON_INLINE #undef TPH_POISSON_MALLOC #undef TPH_POISSON_FREE -/* #undef tph_poisson_vec_invariants #undef tph_poisson_vec_size #undef tph_poisson_vec_append #undef tph_poisson_vec_erase_unordered #undef tph_poisson_vec_shrink_to_fit -*/ #endif// TPH_POISSON_IMPLEMENTATION @@ -1099,15 +1125,11 @@ const tph_poisson_real *tph_poisson_get_samples(tph_poisson_sampling *s) USAGE: - The input points are pruned if - - * There are duplicates points - * The input points are outside of the bounding box (i.e. fail the clipping test function) - * The input points are rejected by the clipper's test function - - The input bounding box is optional (calculated automatically) + Generates a list of samples with the guarantees: (1) No two samples are closer to each other + than some radius; (2) No sample is outside the region bounds. - The input domain is (-FLT_MAX, FLT_MAX] (for floats) + The algorithm tries to fit as many samples as possible into the region without violating the + above requirements. The API consists of these functions: diff --git a/test/src/tph_poisson_d.c b/test/src/tph_poisson_d.c index 2e6e4f93..ee9fa463 100644 --- a/test/src/tph_poisson_d.c +++ b/test/src/tph_poisson_d.c @@ -1,2 +1,7 @@ #define TPH_POISSON_IMPLEMENTATION +/* #define tph_poisson_assert ... */ +/* #define TPH_POISSON_MEMCPY ... */ +/* #define TPH_POISSON_MEMSET ... */ +/* #define TPH_POISSON_MALLOC ... */ +/* #define TPH_POISSON_FREE ... */ #include "tph_poisson_d.h" diff --git a/test/src/tph_poisson_test.hpp b/test/src/tph_poisson_test.hpp deleted file mode 100644 index e69de29b..00000000 diff --git a/thinks/CMakeLists.txt b/thinks/CMakeLists.txt deleted file mode 100644 index c93d8e1c..00000000 --- a/thinks/CMakeLists.txt +++ /dev/null @@ -1,5 +0,0 @@ -# Copyright (C) Tommy Hinks -# This file is subject to the license terms in the LICENSE file -# found in the top-level directory of this distribution. - -add_subdirectory(poisson_disk_sampling) \ No newline at end of file diff --git a/thinks/poisson_disk_sampling/CMakeLists.txt b/thinks/poisson_disk_sampling/CMakeLists.txt deleted file mode 100644 index 96e4b54b..00000000 --- a/thinks/poisson_disk_sampling/CMakeLists.txt +++ /dev/null @@ -1,18 +0,0 @@ -# Copyright (C) Tommy Hinks -# This file is subject to the license terms in the LICENSE file -# found in the top-level directory of this distribution. - -thinks_cc_library( - NAME - poisson_disk_sampling - HDRS - "poisson_disk_sampling.h" - COPTS - ${THINK_DEFAULT_COPTS} - PUBLIC -) - -if (BUILD_TESTING) - add_subdirectory(test) - add_subdirectory(examples) -endif() diff --git a/thinks/poisson_disk_sampling/examples/CMakeLists.txt b/thinks/poisson_disk_sampling/examples/CMakeLists.txt deleted file mode 100644 index 4f7eac41..00000000 --- a/thinks/poisson_disk_sampling/examples/CMakeLists.txt +++ /dev/null @@ -1,62 +0,0 @@ -# Copyright (C) Tommy Hinks -# This file is subject to the license terms in the LICENSE file -# found in the top-level directory of this distribution. - -thinks_cc_executable( - NAME - json_example - SRCS - "json_example.cc" - COPTS - ${THINKS_DEFAULT_COPTS} - DEPS - thinks::poisson_disk_sampling - nlohmann_json::nlohmann_json -) - -thinks_cc_executable( - NAME - periodogram_example - SRCS - "periodogram_example.cc" - COPTS - ${THINKS_DEFAULT_COPTS} - DEPS - thinks::poisson_disk_sampling - d1vanov::simple_fft - nemequ::hedley - nothings::stb -) - -thinks_cc_executable( - NAME - simple_example - SRCS - "simple_example.cc" - COPTS - ${THINKS_DEFAULT_COPTS} - DEPS - thinks::poisson_disk_sampling -) - -thinks_cc_executable( - NAME - vec_traits_in_namespace_example - SRCS - "vec_traits_in_namespace_example.cc" - COPTS - ${THINKS_DEFAULT_COPTS} - DEPS - thinks::poisson_disk_sampling -) - -thinks_cc_executable( - NAME - vec_traits_passed_in_example - SRCS - "vec_traits_passed_in_example.cc" - COPTS - ${THINKS_DEFAULT_COPTS} - DEPS - thinks::poisson_disk_sampling -) diff --git a/thinks/poisson_disk_sampling/examples/json_example.cc b/thinks/poisson_disk_sampling/examples/json_example.cc deleted file mode 100644 index 2b1e7913..00000000 --- a/thinks/poisson_disk_sampling/examples/json_example.cc +++ /dev/null @@ -1,40 +0,0 @@ -// Copyright(C) Tommy Hinks -// This file is subject to the license terms in the LICENSE file -// found in the top-level directory of this distribution. - -#include -#include -#include -#include -#include - -#include "nlohmann/json.hpp" -#include "thinks/poisson_disk_sampling/poisson_disk_sampling.h" - -int main(int /*argc*/, char* /*argv*/[]) { // NOLINT - using json = nlohmann::json; - - constexpr auto kRadius = 3.F; - constexpr std::array kXMin = {-10.F, -10.F}; - constexpr std::array kXMax = {10.F, 10.F}; - constexpr auto kMaxSampleAttempts = std::uint32_t{30}; - constexpr auto kSeed = std::uint32_t{0}; - - const auto samples = thinks::PoissonDiskSampling(kRadius, kXMin, kXMax, - kMaxSampleAttempts, kSeed); - - json j; - j["min"] = kXMin; - j["max"] = kXMax; - j["radius"] = kRadius; - j["samples"] = samples; - - std::ofstream ofs{"./json_example.json"}; - if (!ofs) { - return EXIT_FAILURE; - } - ofs << std::setw(4) << j; - ofs.close(); - - return EXIT_SUCCESS; -} diff --git a/thinks/poisson_disk_sampling/examples/simple_example.cc b/thinks/poisson_disk_sampling/examples/simple_example.cc deleted file mode 100644 index 07593d28..00000000 --- a/thinks/poisson_disk_sampling/examples/simple_example.cc +++ /dev/null @@ -1,30 +0,0 @@ -// Copyright(C) Tommy Hinks -// This file is subject to the license terms in the LICENSE file -// found in the top-level directory of this distribution. - -#include -#include -#include - -#include "thinks/poisson_disk_sampling/poisson_disk_sampling.h" - -int main(int /*argc*/, char* /*argv*/[]) { // NOLINT - constexpr auto kRadius = 2.F; - constexpr auto kXMin = std::array{{-10.F, -10.F}}; - constexpr auto kXMax = std::array{{10.F, 10.F}}; - - // Minimal amount of information provided to sampling function. - const auto samples = thinks::PoissonDiskSampling(kRadius, kXMin, kXMax); - - std::ofstream ofs{"./simple_example.txt"}; - if (!ofs) { - return EXIT_FAILURE; - } - - for (const auto& sample : samples) { - ofs << sample[0] << ", " << sample[1] << '\n'; - } - ofs.close(); - - return EXIT_SUCCESS; -} diff --git a/thinks/poisson_disk_sampling/examples/vec_traits_in_namespace_example.cc b/thinks/poisson_disk_sampling/examples/vec_traits_in_namespace_example.cc deleted file mode 100644 index cdd7463a..00000000 --- a/thinks/poisson_disk_sampling/examples/vec_traits_in_namespace_example.cc +++ /dev/null @@ -1,66 +0,0 @@ -// Copyright(C) Tommy Hinks -// This file is subject to the license terms in the LICENSE file -// found in the top-level directory of this distribution. - -#include -#include -#include - -#include "thinks/poisson_disk_sampling/poisson_disk_sampling.h" - -#if __cplusplus >= 201402L // C++14 or later. -#define PDS_CEXPR constexpr -#else -#define PDS_CEXPR inline -#endif - -namespace { - -struct Vec3 { - float v[3]; // NOLINT -}; - -} // namespace - -namespace thinks { - -// Specialize traits class. -template <> -struct VecTraits { - using ValueType = float; - - static constexpr auto kSize = 3; - - static PDS_CEXPR auto Get(const Vec3& v, const std::size_t i) -> ValueType { - return v.v[i]; - } - - static PDS_CEXPR void Set(Vec3* const v, const std::size_t i, - const ValueType val) { - v->v[i] = val; - } -}; - -} // namespace thinks - -int main(int /*argc*/, char* /*argv*/[]) { // NOLINT - constexpr auto kRadius = 2.F; - constexpr std::array kXMin = {-10.F, -10.F, -10.F}; - constexpr std::array kXMax = {10.F, 10.F, 10.F}; - - // No need to provide traits explicitly. - const auto samples = - thinks::PoissonDiskSampling(kRadius, kXMin, kXMax); - - std::ofstream ofs{"./vec_traits_in_namespace_example.txt"}; - if (!ofs) { - return EXIT_FAILURE; - } - - for (const auto& sample : samples) { - ofs << sample.v[0] << ", " << sample.v[1] << ", " << sample.v[2] << '\n'; - } - ofs.close(); - - return EXIT_SUCCESS; -} diff --git a/thinks/poisson_disk_sampling/examples/vec_traits_passed_in_example.cc b/thinks/poisson_disk_sampling/examples/vec_traits_passed_in_example.cc deleted file mode 100644 index 55258ca2..00000000 --- a/thinks/poisson_disk_sampling/examples/vec_traits_passed_in_example.cc +++ /dev/null @@ -1,61 +0,0 @@ -// Copyright(C) Tommy Hinks -// This file is subject to the license terms in the LICENSE file -// found in the top-level directory of this distribution. - -#include -#include -#include - -#include "thinks/poisson_disk_sampling/poisson_disk_sampling.h" - -#if __cplusplus >= 201402L // C++14 or later. -#define PDS_CEXPR constexpr -#else -#define PDS_CEXPR inline -#endif - -namespace { - -struct Vec3 { - float v[3]; // NOLINT -}; - -// Traits specialized in namespace other than thinks:: -struct Vec3Traits { - using ValueType = float; - - static constexpr auto kSize = 3; - - static PDS_CEXPR auto Get(const Vec3& v, const std::size_t i) -> ValueType { - return v.v[i]; - } - - static PDS_CEXPR void Set(Vec3* const v, const std::size_t i, - const ValueType val) { - v->v[i] = val; - } -}; - -} // namespace - -int main(int /*argc*/, char* /*argv*/[]) { // NOLINT - constexpr auto kRadius = 2.F; - constexpr std::array kXMin = {-10.F, -10.F, -10.F}; - constexpr std::array kXMax = {10.F, 10.F, 10.F}; - - // Traits explicitly provided here. - const auto samples = thinks::PoissonDiskSampling( - kRadius, kXMin, kXMax); - - std::ofstream ofs{"./vec_traits_passed_in_example.txt"}; - if (!ofs) { - return EXIT_FAILURE; - } - - for (const auto& sample : samples) { - ofs << sample.v[0] << ", " << sample.v[1] << ", " << sample.v[2] << '\n'; - } - ofs.close(); - - return EXIT_SUCCESS; -} diff --git a/thinks/poisson_disk_sampling/poisson_disk_sampling.h b/thinks/poisson_disk_sampling/poisson_disk_sampling.h deleted file mode 100644 index 7b835a86..00000000 --- a/thinks/poisson_disk_sampling/poisson_disk_sampling.h +++ /dev/null @@ -1,630 +0,0 @@ -// Copyright(C) Tommy Hinks -// This file is subject to the license terms in the LICENSE file -// found in the top-level directory of this distribution. - -#pragma once - -#include -#include -#include -#include -#include -#include -#include -#include -#include - -#if __cplusplus >= 201402L // C++14 or later. -#define CONSTEXPR14 constexpr -#else -#define CONSTEXPR14 inline -#endif - -namespace thinks { -namespace poisson_disk_sampling_internal { - -// Assumes min_value <= max_value. -template -// NOLINTNEXTLINE -CONSTEXPR14 auto clamped(const ArithT min_value, const ArithT max_value, - const ArithT value) noexcept -> ArithT { - static_assert(std::is_arithmetic::value, "ArithT must be arithmetic"); - return value < min_value ? min_value - : (value > max_value ? max_value : value); -} - -// Returns x squared (not checking for overflow). -template -// NOLINTNEXTLINE -CONSTEXPR14 auto squared(const ArithT x) noexcept -> ArithT { - static_assert(std::is_arithmetic::value, "ArithT must be arithmetic"); - return x * x; -} - -// Returns the squared magnitude of x (not checking for overflow). -template -CONSTEXPR14 auto SquaredMagnitude(const std::array& x) noexcept -> - typename std::array::value_type { - static_assert(std::is_arithmetic::value, "ArithT must be arithmetic"); - constexpr auto kDims = std::tuple_size>::value; - static_assert(kDims >= 1, "dimensions must be >= 1"); - - auto m = squared(x[0]); - for (std::size_t i = 1; i < kDims; ++i) { - m += squared(x[i]); - } - return m; -} - -// Returns the squared distance between the vectors u and v. -template -CONSTEXPR14 auto SquaredDistance(const VecT& u, const VecT& v) noexcept -> - typename VecTraitsT::ValueType { - static_assert(VecTraitsT::kSize >= 1, "dimensions must be >= 1"); - - auto d = squared(VecTraitsT::Get(u, 0) - VecTraitsT::Get(v, 0)); - for (std::size_t i = 1; i < VecTraitsT::kSize; ++i) { - d += squared(VecTraitsT::Get(u, i) - VecTraitsT::Get(v, i)); - } - return d; -} - -// Returns true if x is element-wise inclusively inside x_min and x_max, -// otherwise false. Assumes that x_min is element-wise less than x_max. -template -CONSTEXPR14 auto InsideBounds(const VecT& sample, - const std::array& x_min, - const std::array& x_max) noexcept - -> bool { - constexpr auto kDims = std::tuple_size>::value; - static_assert(VecTraitsT::kSize == kDims, "dimensionality mismatch"); - - for (std::size_t i = 0; i < VecTraitsT::kSize; ++i) { - const auto xi = static_cast(VecTraitsT::Get(sample, i)); - if (!(x_min[i] <= xi && xi <= x_max[i])) { - return false; - } - } - return true; -} - -// Erase the value at given index in the vector. The vector is -// guaranteed to decrease in size by one. Note that the ordering of elements -// may change as a result of calling this function. -// -// Assumes that the vector argument is non-null and that the index is valid. -// (Cannot be called on an empty vector since no valid indices exist) -template -void EraseUnordered(std::vector* const v, const std::size_t index) noexcept { - // Element at index gets same value as last element. - // Last element is removed, resulting in a vector that - // has one fewer elements with the value that was previously - // at index. - (*v)[index] = v->back(); - v->pop_back(); // O(1). -} - -// Increment index such that repeated calls to this function enumerate -// all iterations between min_index and max_index (inclusive). -// Assumes that min_index is element-wise less than or equal to max_index. -template -CONSTEXPR14 auto Iterate(const std::array& min_index, - const std::array& max_index, - std::array* const index) noexcept -> bool { - static_assert(std::is_integral::value, "IntT must be integral"); - constexpr auto kDims = std::tuple_size>::value; - static_assert(kDims >= 1, "dimensions must be >= 1"); - - std::size_t i = 0; - for (; i < kDims; ++i) { - (*index)[i]++; - if ((*index)[i] <= max_index[i]) { - break; - } - (*index)[i] = min_index[i]; - } - return i == kDims ? false : true; -} - -// Stateless and repeatable function that returns a -// pseduo-random number in the range [0, 0xFFFFFFFF]. -CONSTEXPR14 auto Hash(const std::uint32_t seed) noexcept -> std::uint32_t { - // So that we can use unsigned int literals, e.g. 42u. - static_assert(sizeof(unsigned int) == sizeof(std::uint32_t), - "integer size mismatch"); - - auto i = std::uint32_t{(seed ^ 12345391U) * 2654435769U}; // NOLINT - i ^= (i << 6U) ^ (i >> 26U); // NOLINT - i *= 2654435769U; // NOLINT - i += (i << 5U) ^ (i >> 12U); // NOLINT - return i; -} - -// Returns a pseduo-random number in the range [0, 0xFFFFFFFF]. -// Note that seed is incremented for each invokation. -CONSTEXPR14 auto Rand(std::uint32_t* const seed) noexcept -> std::uint32_t { - // Not worrying about seed "overflow" since it is unsigned. - return Hash((*seed)++); -} - -// Returns a pseduo-random number in the range [0, 1]. -template -CONSTEXPR14 auto NormRand(std::uint32_t* const seed) noexcept -> FloatT { - static_assert(std::is_floating_point::value, - "FloatT must be floating point"); - // TODO(thinks): clamped? - return (1 / static_cast(std::numeric_limits::max())) * - static_cast(Rand(seed)); -} - -// Returns a pseduo-random number in the range [offset, offset + range]. -// Assumes range > 0. -template -CONSTEXPR14 auto RangeRand(const FloatT offset, const FloatT range, - std::uint32_t* const seed) noexcept -> FloatT { - return offset + range * NormRand(seed); -} - -// Returns an array where each element has been assigned using RangeRand, -// with bounds taken from the corresponding element in x_min and x_max. -// Assumes that x_min[i] < x_max[i]. -template -CONSTEXPR14 auto ArrayRangeRand(const std::array& x_min, - const std::array& x_max, - std::uint32_t* const seed) noexcept - -> std::array { - constexpr auto kDims = std::tuple_size>::value; - - std::array a = {}; - for (std::size_t i = 0; i < kDims; ++i) { - a[i] = - RangeRand(/* offset */ x_min[i], /* range */ x_max[i] - x_min[i], seed); - } - return a; -} - -// Returns an array where each element has been assigned using RangeRand, -// with bounds taken from x_min and x_max. -// Assumes that x_min < x_max. -template -CONSTEXPR14 auto ArrayRangeRand(const FloatT x_min, const FloatT x_max, - std::uint32_t* const seed) noexcept - -> std::array { - constexpr auto kDims = std::tuple_size>::value; - - const auto range = x_max - x_min; - std::array a = {}; - for (std::size_t i = 0; i < kDims; ++i) { - a[i] = RangeRand(/* offset */ x_min, range, seed); - } - return a; -} - -// See ArrayRangeRand. -template -CONSTEXPR14 auto VecRangeRand(const std::array& x_min, - const std::array& x_max, - std::uint32_t* const seed) noexcept -> VecT { - constexpr auto kDims = std::tuple_size>::value; - static_assert(VecTraitsT::kSize == kDims, "dimensionality mismatch"); - - VecT v = {}; - const auto a = ArrayRangeRand(x_min, x_max, seed); - for (std::size_t i = 0; i < kDims; ++i) { - VecTraitsT::Set(&v, i, static_cast(a[i])); - } - return v; -} - -// Returns a pseudo-random index in the range [0, size - 1]. -CONSTEXPR14 auto IndexRand(const std::size_t size, - std::uint32_t* const seed) noexcept -> std::size_t { - constexpr auto kEps = 0.0001F; - return static_cast( - RangeRand(float{0}, static_cast(size) - kEps, seed)); -} - -template -class Grid { - public: - using CellType = std::int32_t; - using IndexType = std::array; - - static constexpr auto kDims = std::tuple_size::value; - - static_assert(kDims >= 1, "grid dimensionality must be >= 1"); - static_assert(std::is_floating_point::value, - "FloatT must be floating point"); - - explicit Grid(const FloatT sample_radius, const std::array& x_min, - const std::array& x_max) noexcept - : sample_radius_(sample_radius), - dx_(GetDx_(sample_radius_)), - dx_inv_(FloatT{1} / dx_), - x_min_(x_min), - x_max_(x_max), - size_(GetGridSize_(x_min_, x_max_, dx_inv_)), - cells_(GetCells_(size_)) {} - - auto sample_radius() const noexcept -> FloatT { // NOLINT - return sample_radius_; - } - - auto size() const noexcept -> IndexType { return size_; } // NOLINT - - // Returns the index for a position along the i'th axis. - // Note that the returned index may be negative. - template - // NOLINTNEXTLINE - auto AxisIndex(const std::size_t i, const FloatT2 pos) const noexcept -> - typename IndexType::value_type { - using IndexValueType = typename IndexType::value_type; - - return static_cast((static_cast(pos) - x_min_[i]) * - dx_inv_); - } - - // Note that the returned index elements may be negative. - template - // NOLINTNEXTLINE - auto IndexFromSample(const VecT& sample) const noexcept -> IndexType { - static_assert(VecTraitsT::kSize == kDims, "dimensionality mismatch"); - - IndexType index = {}; - for (std::size_t i = 0; i < kDims; ++i) { - index[i] = AxisIndex(i, VecTraitsT::Get(sample, i)); - } - return index; - } - - // NOLINTNEXTLINE - auto Cell(const IndexType& index) const noexcept -> CellType { - return cells_[LinearIndex_(index)]; - } - - // NOLINTNEXTLINE - auto Cell(const IndexType& index) noexcept -> CellType& { - return cells_[LinearIndex_(index)]; - } - - private: - FloatT sample_radius_; - FloatT dx_; - FloatT dx_inv_; - std::array x_min_; - std::array x_max_; - IndexType size_; - std::vector cells_; - - // Assumes that all elements in index are >= 0. - // NOLINTNEXTLINE - auto LinearIndex_(const IndexType& index) const noexcept -> std::size_t { - auto k = static_cast(index[0]); - auto d = std::size_t{1}; - for (auto i = std::size_t{1}; i < kDims; ++i) { - // Note: Not checking for "overflow". - d *= static_cast(size_[i - 1]); - k += static_cast(index[i]) * d; - } - return k; - } - - // Assumes sample_radius is > 0. - static auto GetDx_(const FloatT sample_radius) noexcept -> FloatT { - // The grid cell size should be such that each cell can only - // contain one sample. We apply a scaling factor to avoid - // numerical issues. - constexpr auto kEps = static_cast(0.001); - constexpr auto kScale = FloatT{1} - kEps; - return kScale * sample_radius / std::sqrt(static_cast(N)); - } - - // Assumes that x_min is element-wise less than x_max. - // Assumes dx_inv > 0. - static auto GetGridSize_(const std::array& x_min, - const std::array& x_max, - const FloatT dx_inv) noexcept -> IndexType { - // Compute size in each dimension using grid cell size (dx). - auto grid_size = IndexType{}; - for (auto i = std::size_t{0}; i < kDims; ++i) { - grid_size[i] = static_cast( - std::ceil((x_max[i] - x_min[i]) * dx_inv)); - } - return grid_size; - } - - static auto GetCells_(const IndexType& size) noexcept - -> std::vector { - // Initialize cells with value -1, indicating no sample there. - // Cell values are later set to indices of samples. - const auto linear_size = - std::accumulate(std::begin(size), std::end(size), std::size_t{1}, - std::multiplies()); - return std::vector(linear_size, -1); - } -}; - -// Named contructor to help with type deduction. -template -auto MakeGrid(const FloatT sample_radius, const std::array& x_min, - const std::array& x_max) noexcept -> Grid { - return Grid(sample_radius, x_min, x_max); -} - -template -CONSTEXPR14 auto ValidBounds(const std::array& x_min, - const std::array& x_max) noexcept - -> bool { - static_assert(std::is_arithmetic::value, "type must be arithmetic"); - constexpr auto kDims = std::tuple_size>::value; - static_assert(kDims >= 1, "bounds dimensionality must be >= 1"); - - for (std::size_t i = 0; i < kDims; ++i) { - if (!(x_max[i] > x_min[i])) { - return false; - } - } - return true; -} - -template -struct ActiveSample { - VecT position; - std::size_t active_index; - std::uint32_t sample_index; -}; - -// Returns a pseudo-randomly selected active sample from the active indices. -template -auto RandActiveSample(const std::vector& active_indices, - const std::vector& samples, - std::uint32_t* const seed) noexcept - -> ActiveSample { - ActiveSample active_sample = {}; - active_sample.active_index = IndexRand(active_indices.size(), seed); - active_sample.sample_index = active_indices[active_sample.active_index]; - active_sample.position = samples[active_sample.sample_index]; - return active_sample; -} - -// Returns a pseudo-random coordinate that is guaranteed be at a -// distance [radius, 2 * radius] from center. -template -CONSTEXPR14 auto RandAnnulusSample(const VecT& center, const FloatT radius, - std::uint32_t* const seed) noexcept -> VecT { - VecT p = {}; - for (;;) { - // Generate a random component in the range [-2, 2] for each dimension. - const auto offset = - ArrayRangeRand(FloatT{-2}, FloatT{2}, seed); - - // The randomized offset is not guaranteed to be within the radial - // distance that we need to guarantee. If we found an offset with - // magnitude in the range (1, 2] we are done, otherwise generate a new - // offset. Continue until a valid offset is found. - const auto r2 = SquaredMagnitude(offset); - if (FloatT{1} < r2 && r2 <= FloatT{4}) { - // Found a valid offset. - // Add the offset scaled by radius to the center coordinate to - // produce the final sample. - for (std::size_t i = 0; i < VecTraitsT::kSize; ++i) { - VecTraitsT::Set(&p, i, - static_cast( - static_cast(VecTraitsT::Get(center, i)) + - radius * offset[i])); - } - break; - } - } - return p; -} - -// Add a new sample and update data structures. -template -void AddSample(const VecT& sample, std::vector* const samples, - std::vector* const active_indices, - Grid* const grid) noexcept { - const auto sample_index = samples->size(); // New sample index. - samples->push_back(sample); - active_indices->push_back(static_cast(sample_index)); - const auto grid_index = grid->template IndexFromSample(sample); - grid->Cell(grid_index) = - static_cast(sample_index); // No range check! -} - -template -struct GridIndexRange { - IndexT min_index; - IndexT max_index; -}; - -// Returns the grid neighborhood of the sample using the radius of the grid. -template -auto GridNeighborhood(const VecT& sample, const Grid& grid) noexcept - -> GridIndexRange::IndexType> { - using GridIndexType = typename std::decay::type::IndexType; - using GridIndexValueType = typename GridIndexType::value_type; - - constexpr auto kDims = std::tuple_size::value; - static_assert(VecTraitsT::kSize == kDims, "dimensionality mismatch"); - - auto min_index = GridIndexType{}; - auto max_index = GridIndexType{}; - const auto grid_size = grid.size(); - const auto radius = grid.sample_radius(); - for (auto i = std::size_t{0}; i < kDims; ++i) { - const auto xi_min = GridIndexValueType{0}; - const auto xi_max = static_cast(grid_size[i] - 1); - const auto xi = static_cast(VecTraitsT::Get(sample, i)); - const auto xi_sub = grid.AxisIndex(i, xi - radius); - const auto xi_add = grid.AxisIndex(i, xi + radius); - min_index[i] = clamped(xi_min, xi_max, xi_sub); - max_index[i] = clamped(xi_min, xi_max, xi_add); - } - return GridIndexRange{min_index, max_index}; -} - -// Returns true if there exists another sample within the radius used to -// construct the grid, otherwise false. -template -auto ExistingSampleWithinRadius( - const VecT& sample, const std::uint32_t active_sample_index, - const std::vector& samples, const Grid& grid, - const typename Grid::IndexType& min_index, - const typename Grid::IndexType& max_index) noexcept -> bool { - auto index = min_index; - const auto r_squared = squared(grid.sample_radius()); - do { - const auto cell_index = grid.Cell(index); - if (cell_index >= 0 && - static_cast(cell_index) != active_sample_index) { - const auto cell_sample = samples[static_cast(cell_index)]; - const auto d = - static_cast(SquaredDistance(sample, cell_sample)); - if (d < r_squared) { - return true; - } - } - } while (Iterate(min_index, max_index, &index)); - - return false; -} - -} // namespace poisson_disk_sampling_internal - -// Generic template for vector traits. Users may specialize this template -// for their own classes. -// -// Specializations must have the following static interface, here using -// an example type 'MyVec'. -// -// struct MyVecTraits -// { -// using ValueType = ... -// static constexpr std::size_t kSize = ... -// -// static ValueType Get(const MyVec&, std::size_t) -// static void Get(MyVec* const, std::size_t, ValueType) -// } -// -// See specialization for std::array below as an example. -template -struct VecTraits; - -// Specialization of vector traits for std::array. -template -struct VecTraits> { - using ValueType = typename std::array::value_type; - static_assert(std::is_floating_point::value, - "ValueType must be floating point"); - - static constexpr auto kSize = std::tuple_size>::value; - static_assert(kSize >= 1, "kSize must be >= 1"); - - // No bounds checking. - static CONSTEXPR14 auto Get(const std::array& vec, - const std::size_t i) noexcept -> ValueType { - return vec[i]; - } - - // No bounds checking. - static CONSTEXPR14 void Set(std::array* const vec, - const std::size_t i, - const ValueType value) noexcept { - (*vec)[i] = value; - } -}; - -// Returns a list of samples with the guarantees: -// - No two samples are closer to each other than radius. -// - No sample is outside the region [x_min, x_max]. -// -// The algorithm tries to fit as many samples as possible -// into the region without violating the above requirements. -// -// If the arguments are invalid an empty vector is returned. -// The arguments are invalid if: -// - Radius is <= 0, or -// - x_min[i] >= x_max[i], or -// - max_sample_attempts == 0. -template , - typename VecTraitsT = VecTraits> -auto PoissonDiskSampling(const FloatT radius, - const std::array& x_min, - const std::array& x_max, - const std::uint32_t max_sample_attempts = 30, - const std::uint32_t seed = 0) noexcept - -> std::vector { - namespace pds = poisson_disk_sampling_internal; - - // Validate input. - if (!(radius > FloatT{0}) || !(max_sample_attempts > 0) || - !pds::ValidBounds(x_min, x_max)) { - // Returning an empty list of samples indicates an error, - // since for any valid input there is always at least one sample. - return std::vector{}; - } - - // Acceleration grid. - auto grid = pds::MakeGrid(radius, x_min, x_max); - - auto samples = std::vector{}; - auto active_indices = std::vector{}; - auto local_seed = seed; - - // Add first sample randomly within bounds. - // No need to check (non-existing) neighbors. - pds::AddSample( - pds::VecRangeRand(x_min, x_max, &local_seed), &samples, - &active_indices, &grid); - - while (!active_indices.empty()) { - // Randomly choose an active sample. A sample is considered active - // until failed attempts have been made to generate a new sample within - // its annulus. - const auto active_sample = - pds::RandActiveSample(active_indices, samples, &local_seed); - - auto attempt_count = std::uint32_t{0}; - while (attempt_count < max_sample_attempts) { - // Randomly create a candidate sample inside the active sample's annulus. - const auto cand_sample = pds::RandAnnulusSample( - active_sample.position, grid.sample_radius(), &local_seed); - - // Check if candidate sample is within bounds. - if (pds::InsideBounds(cand_sample, x_min, x_max)) { - // Check candidate sample proximity to nearby samples. - const auto grid_neighbors = - pds::GridNeighborhood(cand_sample, grid); - const auto existing_sample = - pds::ExistingSampleWithinRadius( - cand_sample, active_sample.sample_index, samples, grid, - grid_neighbors.min_index, grid_neighbors.max_index); - if (!existing_sample) { - // No existing samples where found to be too close to the - // candidate sample, no further attempts necessary. - pds::AddSample(cand_sample, &samples, &active_indices, - &grid); - break; - } - // Else: The candidate sample is too close to an existing sample, - // move on to next attempt. - } - // Else: The candidate sample is out-of-bounds. - - ++attempt_count; - } - - if (attempt_count == max_sample_attempts) { - // No valid sample was found on the disk of the active sample, - // remove it from the active list. - pds::EraseUnordered(&active_indices, active_sample.active_index); - } - } - - return samples; -} - -} // namespace thinks - -#undef CONSTEXPR14 diff --git a/thinks/poisson_disk_sampling/test/CMakeLists.txt b/thinks/poisson_disk_sampling/test/CMakeLists.txt deleted file mode 100644 index 93e6cac2..00000000 --- a/thinks/poisson_disk_sampling/test/CMakeLists.txt +++ /dev/null @@ -1,17 +0,0 @@ -# Copyright (C) Tommy Hinks -# This file is subject to the license terms in the LICENSE file -# found in the top-level directory of this distribution. - -thinks_cc_test( - NAME - poisson_disk_sampling_test - SRCS - "catch_main.cc" - "poisson_disk_sampling_test.cc" - COPTS - ${THINKS_TEST_COPTS} - DEPS - thinks::poisson_disk_sampling - Catch2::Catch2 - ${CMAKE_THREAD_LIBS_INIT} -) diff --git a/thinks/poisson_disk_sampling/test/catch_main.cc b/thinks/poisson_disk_sampling/test/catch_main.cc deleted file mode 100644 index d89af05f..00000000 --- a/thinks/poisson_disk_sampling/test/catch_main.cc +++ /dev/null @@ -1,6 +0,0 @@ -// Copyright(C) Tommy Hinks -// This file is subject to the license terms in the LICENSE file -// found in the top-level directory of this distribution. - -#define CATCH_CONFIG_MAIN -#include "catch2/catch.hpp" diff --git a/thinks/poisson_disk_sampling/test/poisson_disk_sampling_test.cc b/thinks/poisson_disk_sampling/test/poisson_disk_sampling_test.cc deleted file mode 100644 index 1210d76a..00000000 --- a/thinks/poisson_disk_sampling/test/poisson_disk_sampling_test.cc +++ /dev/null @@ -1,449 +0,0 @@ -// Copyright(C) Tommy Hinks -// This file is subject to the license terms in the LICENSE file -// found in the top-level directory of this distribution. - -#include "thinks/poisson_disk_sampling/poisson_disk_sampling.h" - -#include -#include -#include -#include -#include -#include -#include - -#include "catch2/catch.hpp" - -#if __cplusplus >= 201402L // C++14 or later. -#define CONSTEXPR14 constexpr -#else -#define CONSTEXPR14 inline -#endif - -namespace { - -// Simple small vector type. -template -struct Vec { - std::array m; -}; - -} // namespace - -namespace thinks { - -// Specializations of VecTraits in thinks:: namespace. -template -struct VecTraits> { - using ValueType = typename decltype(Vec::m)::value_type; - - static_assert(sizeof(Vec) == N * sizeof(T), - "Vec type must be tightly packed"); - static constexpr auto kSize = sizeof(Vec) / sizeof(ValueType); - - static constexpr auto Get(const Vec& v, const std::size_t i) noexcept - -> ValueType { - return v.m[i]; - } - - static CONSTEXPR14 void Set(Vec* const v, const std::size_t i, - const ValueType val) noexcept { - v->m[i] = val; - } -}; - -namespace { - -auto ThreadCount(const std::size_t max_thread_count) noexcept -> std::size_t { - return std::thread::hardware_concurrency() > 0 - ? std::min(static_cast( - std::thread::hardware_concurrency()), - max_thread_count) - : 1; -} - -template -constexpr auto squared(const T x) noexcept -> T { // NOLINT - return x * x; -} - -// Returns the squared distance between u and v. Not checking for overflow. -template -CONSTEXPR14 auto SquaredDistance(const VecT& u, const VecT& v) noexcept -> - typename VecTraitsT::ValueType { - static_assert(VecTraitsT::kSize >= 1, "vec dimensionality must be >= 1"); - - auto d = squared(VecTraitsT::Get(u, 0) - VecTraitsT::Get(v, 0)); - for (std::size_t i = 1; i < VecTraitsT::kSize; ++i) { - d += squared(VecTraitsT::Get(u, i) - VecTraitsT::Get(v, i)); - } - return d; -} - -// Brute-force (with some tricks) verification that the distance between each -// possible sample pair meets the Poisson requirement, i.e. is greater than some -// radius. Returns true if Poisson criteria is met for all samples, otherwise -// false. -template -auto VerifyPoisson(const std::vector& samples, - const FloatT radius) noexcept -> bool { - if (samples.empty()) { - return false; - } - if (samples.size() == 1) { - return true; - } - - // Setup threading. - // Avoid spawning more threads than there are samples (although very - // unlikely). - const auto sample_count = samples.size(); - const auto thread_count = ThreadCount(sample_count); - - // Launch threads. - std::vector> futures; - for (std::size_t i = 0; i < thread_count; ++i) { - // Each thread works on local copy of samples. - futures.emplace_back(std::async(std::launch::async, [=]() { - // We know that distance is symmetrical, such that - // dist(s[j], s[k]) == dist(s[k], s[j]). Therefore - // we need only compute the upper half of the matrix (excluding diagonal). - // - // Load balance threads such that "short" (small j) and "long" (large j) - // columns are divided evenly among threads. - const auto r_sqr = squared(radius); - for (std::size_t j = i; j < sample_count; j += thread_count) { - const auto sj = samples[j]; - const auto k_max = j; - for (std::size_t k = 0; k < k_max; ++k) { - const auto sk = samples[k]; - const auto dist_sqr = static_cast( - SquaredDistance>(sj, sk)); - - // Fail for NaN. - if (!(dist_sqr > r_sqr)) { - return false; - } - } - } - return true; - })); - } - - // clang-format off - // Check results. - for (auto&& f : futures) { f.wait(); } - return std::all_of(std::begin(futures), std::end(futures), - [](std::future& f) { return f.get(); }); - // clang-format on -} - -template -CONSTEXPR14 auto SampleInsideBounds(const VecT& sample, - const std::array& x_min, - const std::array& x_max) noexcept - -> bool { - using VecTraitsType = thinks::VecTraits; - - constexpr auto kDims = std::tuple_size>::value; - static_assert(kDims == VecTraitsType::kSize, - "sample/bounds dimensionality mismatch"); - - for (std::size_t i = 0; i < kDims; ++i) { - const auto xi = static_cast(VecTraitsType::Get(sample, i)); - - // Fail if x_min > x_max. - if (!(x_min[i] <= xi && xi <= x_max[i])) { - return false; - } - } - return true; -} - -// Returns true if all samples are within the bounds specified by x_min and -// x_max, otherwise false. -template -auto VerifyBounds(const std::vector& samples, - const std::array& x_min, - const std::array& x_max) noexcept -> bool { - if (samples.empty()) { - return false; - } - - // Setup threading. - // Avoid spawning more threads than there are samples (although very - // unlikely). - const auto sample_count = samples.size(); - const auto thread_count = ThreadCount(sample_count); - const auto chunk_size = - sample_count / thread_count + (sample_count % thread_count != 0); - - // Launch threads. - std::vector> futures; - for (std::size_t i = 0; i < thread_count; ++i) { - // Each thread works on local copy of samples. - futures.emplace_back(std::async(std::launch::async, [=]() { - const std::size_t begin = i * chunk_size; - const std::size_t end = std::min((i + 1) * chunk_size, sample_count); - for (std::size_t j = begin; j < end; ++j) { - if (!SampleInsideBounds(samples[j], x_min, x_max)) { - return false; - } - } - return true; - })); - } - - // clang-format off - // Check results. - for (auto&& f : futures) { f.wait(); } - return std::all_of(std::begin(futures), std::end(futures), - [](std::future& f) { return f.get(); }); - // clang-format on -} - -template > -auto VerifyPoissonDiskSampling(const FloatT radius, - const std::array& x_min, - const std::array& x_max, - const std::uint32_t max_sample_attempts = 30, - const std::uint32_t seed = 0) noexcept -> bool { - FloatT config_radius = radius; - -#ifndef NDEBUG - // Reduce the number of samples for debug builds to decrease - // verification times. Larger radius gives fewer samples. - config_radius *= 2; -#endif - - const auto samples = thinks::PoissonDiskSampling( - config_radius, x_min, x_max, max_sample_attempts, seed); - return VerifyBounds(samples, x_min, x_max) && VerifyPoisson(samples, radius); -} - -template -constexpr auto SampleTestRadius() noexcept -> FloatT { - return FloatT{2}; -} - -template -struct SampleTestBounds; - -// clang-format off -template -struct SampleTestBounds { - static constexpr auto x_min() noexcept -> std::array { return {{-100, -100}}; } // NOLINT - static constexpr auto x_max() noexcept -> std::array { return {{100, 100}}; } // NOLINT -}; - -template -struct SampleTestBounds { - static constexpr auto x_min() noexcept -> std::array { return {{-20, -20, -20}}; } // NOLINT - static constexpr auto x_max() noexcept -> std::array { return {{20, 20, 20}}; } // NOLINT -}; - -template -struct SampleTestBounds { - static constexpr auto x_min() noexcept -> std::array { return {{-10, -10, -10, -10}}; } // NOLINT - static constexpr auto x_max() noexcept -> std::array { return {{10, 10, 10, 10}}; } // NOLINT -}; -// clang-format on - -TEST_CASE("Verify samples ", "[container]") { - SECTION("N = 2") { - REQUIRE(VerifyPoissonDiskSampling(SampleTestRadius(), - SampleTestBounds::x_min(), - SampleTestBounds::x_max())); - REQUIRE(VerifyPoissonDiskSampling(SampleTestRadius(), - SampleTestBounds::x_min(), - SampleTestBounds::x_max())); - } - SECTION("N = 3") { - REQUIRE(VerifyPoissonDiskSampling(SampleTestRadius(), - SampleTestBounds::x_min(), - SampleTestBounds::x_max())); - REQUIRE(VerifyPoissonDiskSampling(SampleTestRadius(), - SampleTestBounds::x_min(), - SampleTestBounds::x_max())); - } - SECTION("N = 4") { - REQUIRE(VerifyPoissonDiskSampling(SampleTestRadius(), - SampleTestBounds::x_min(), - SampleTestBounds::x_max())); - REQUIRE(VerifyPoissonDiskSampling(SampleTestRadius(), - SampleTestBounds::x_min(), - SampleTestBounds::x_max())); - } -} - -TEST_CASE("Verify samples ", "[container]") { - SECTION("N = 2") { - using VecType = Vec; - REQUIRE(VerifyPoissonDiskSampling( - SampleTestRadius(), SampleTestBounds::x_min(), - SampleTestBounds::x_max())); - REQUIRE(VerifyPoissonDiskSampling( - SampleTestRadius(), SampleTestBounds::x_min(), - SampleTestBounds::x_max())); - } - SECTION("N = 3") { - using VecType = Vec; - REQUIRE(VerifyPoissonDiskSampling( - SampleTestRadius(), SampleTestBounds::x_min(), - SampleTestBounds::x_max())); - REQUIRE(VerifyPoissonDiskSampling( - SampleTestRadius(), SampleTestBounds::x_min(), - SampleTestBounds::x_max())); - } - SECTION("N = 4") { - using VecType = Vec; - REQUIRE(VerifyPoissonDiskSampling( - SampleTestRadius(), SampleTestBounds::x_min(), - SampleTestBounds::x_max())); - REQUIRE(VerifyPoissonDiskSampling( - SampleTestRadius(), SampleTestBounds::x_min(), - SampleTestBounds::x_max())); - } -} - -TEST_CASE("Verify max sampling attempts") { - // Verify that we get a denser sampling, i.e. more samples, - // when we increase the max sample attempts parameter (with - // all other parameters constant). - - constexpr auto kRadius = 0.5F; - constexpr auto kXMin = std::array{{-10.F, -10.F}}; - constexpr auto kXMax = std::array{{10.F, 10.F}}; - - const auto samples_10 = thinks::PoissonDiskSampling( - kRadius, kXMin, kXMax, /* max_sample_attempts */ 10); - const auto samples_40 = thinks::PoissonDiskSampling( - kRadius, kXMin, kXMax, /* max_sample_attempts */ 40); - - REQUIRE(samples_10.size() < samples_40.size()); -} - -TEST_CASE("Verify seed") { - // Verify that different seeds give different sample distributions. - - constexpr auto kRadius = 0.5F; - constexpr auto kXMin = std::array{{-10.F, -10.F}}; - constexpr auto kXMax = std::array{{10.F, 10.F}}; - constexpr auto kMaxSampleAttempts = std::uint32_t{20}; - - const auto samples_1981 = thinks::PoissonDiskSampling( - kRadius, kXMin, kXMax, kMaxSampleAttempts, /* seed */ 1981); - const auto samples_1337 = thinks::PoissonDiskSampling( - kRadius, kXMin, kXMax, kMaxSampleAttempts, /* seed */ 1337); - - // For each sample in the first point set compute the smallest - // distance checking every sample in the second point set. - // Then, if the smallest distance is larger than some threshold - // we say that the sample from the first point set is distinct - // from every sample in the second point set. Thus the two - // distributions must be different. - const auto sample_count_1981 = samples_1981.size(); - const auto sample_count_1337 = samples_1337.size(); - auto distinct_sample_found = false; - for (std::size_t i = 0; i < sample_count_1981; ++i) { - const auto si = samples_1981[i]; - auto min_sqr_dist = std::numeric_limits::max(); - for (std::size_t j = 0; j < sample_count_1337; ++j) { - const auto sj = samples_1981[j]; - min_sqr_dist = std::min( - min_sqr_dist, - SquaredDistance>>(si, sj)); - } - if (min_sqr_dist > 0.1F) { // NOLINT - distinct_sample_found = true; - break; - } - } - - REQUIRE(distinct_sample_found); -} - -struct ValidBounds { - static constexpr auto XMin() noexcept -> std::array { - return {{-1, -1}}; - } - static constexpr auto XMax() noexcept -> std::array { - return {{1, 1}}; - } -}; - -TEST_CASE("Invalid arguments", "[container]") { - constexpr auto kValidRadius = 1.F; - - SECTION("radius == 0") { - const auto samples = thinks::PoissonDiskSampling( - /* radius */ 0.F, // NOLINT - ValidBounds::XMin(), ValidBounds::XMax()); - REQUIRE(samples.empty()); - } - - SECTION("radius < 0") { - const auto samples = thinks::PoissonDiskSampling( - /* radius */ -1.F, // NOLINT - ValidBounds::XMin(), ValidBounds::XMax()); - REQUIRE(samples.empty()); - } - - SECTION("x_min == x_max") { - { - const auto samples = thinks::PoissonDiskSampling( - kValidRadius, - /* x_min */ std::array{{1.F, 1.F}}, // NOLINT - /* x_max */ std::array{{1.F, 1.F}}); // NOLINT - REQUIRE(samples.empty()); - } - { - const auto samples = thinks::PoissonDiskSampling( - kValidRadius, - /* x_min */ std::array{{1.F, -1.F}}, // NOLINT - /* x_max */ std::array{{1.F, 1.F}}); // NOLINT - REQUIRE(samples.empty()); - } - { - const auto samples = thinks::PoissonDiskSampling( - kValidRadius, - /* x_min */ std::array{{-1.F, 1.F}}, // NOLINT - /* x_max */ std::array{{1.F, 1.F}}); // NOLINT - REQUIRE(samples.empty()); - } - } - - SECTION("x_min > x_max") { - { - const auto samples = thinks::PoissonDiskSampling( - kValidRadius, - /* x_min */ std::array{{1.F, 1.F}}, // NOLINT - /* x_max */ std::array{{-1.F, -1.F}}); // NOLINT - REQUIRE(samples.empty()); - } - { - const auto samples = thinks::PoissonDiskSampling( - kValidRadius, - /* x_min */ std::array{{1.F, -1.F}}, // NOLINT - /* x_max */ std::array{{-1.F, 1.F}}); // NOLINT - REQUIRE(samples.empty()); - } - { - const auto samples = thinks::PoissonDiskSampling( - kValidRadius, - /* x_min */ std::array{{-1.F, 1.F}}, // NOLINT - /* x_max */ std::array{{1.F, -1.F}}); // NOLINT - REQUIRE(samples.empty()); - } - } - - SECTION("max_sample_attempts == 0") { - const auto samples = thinks::PoissonDiskSampling( - kValidRadius, ValidBounds::XMin(), ValidBounds::XMax(), - /* max_sample_attempts */ std::uint32_t{0}); // NOLINT - REQUIRE(samples.empty()); - } -} - -} // namespace -} // namespace thinks