nerdcruft/libcruft-util
1
0
Fork 0
Code Issues Releases Activity

geom/sample: align the interfaces for volume/surface sampling

Browse Source
This commit is contained in:
Danny Robson 2018-11-26 14:08:50 +11:00
parent 38906862a3
commit 1e0a87d5df
11 changed files with 315 additions and 206 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

6
CMakeLists.txt
View File

@ -276,6 +276,7 @@ list (
geom/iostream.hpp
geom/line.hpp
geom/line.cpp
geom/ops.cpp
geom/ops.hpp
geom/plane.cpp
geom/plane.hpp
@ -283,8 +284,9 @@ list (
geom/ray.hpp
geom/rect.cpp
geom/rect.hpp
geom/sample.cpp
geom/sample.hpp
geom/sample/fwd.hpp
geom/sample/surface.hpp
geom/sample/volume.hpp
geom/segment.cpp
geom/segment.hpp
geom/sphere.cpp

24
geom/aabb.hpp
View File

@ -138,23 +138,35 @@ namespace cruft::geom {
///////////////////////////////////////////////////////////////////////////////
#include "./sample.hpp"
#include "sample/fwd.hpp"
#include <random>
namespace cruft::geom {
namespace cruft::geom::sample {
template <size_t S, typename T>
struct sampler<aabb<S,T>> {
class volume<aabb<S,T>> {
public:
using shape_type = aabb<S,T>;
volume (shape_type const &target):
m_target (target)
{ ; }
template <typename GeneratorT>
static auto
eval (aabb<S,T> shape, GeneratorT &g)
auto
eval (GeneratorT &&g) const noexcept
{
point<S,T> res;
for (size_t i = 0; i < S; ++i)
res[i] = cruft::random::uniform (shape.lo[i], shape.hi[i], g);
res[i] = cruft::random::uniform (m_target.lo[i], m_target.hi[i], g);
return res;
}
private:
shape_type const &m_target;
};
}

29
geom/ellipse.hpp
View File

@ -138,20 +138,29 @@ namespace cruft::geom {
///////////////////////////////////////////////////////////////////////////////
#include "sample.hpp"
#include "sample/fwd.hpp"
#include "../random.hpp"
#include <cmath>
#include <random>
namespace cruft::geom {
namespace cruft::geom::sample {
/// Specialisation for uniform sampling of ellipses
template <typename T>
struct sampler<ellipse<2,T>>
{
class volume<ellipse<2,T>> {
public:
using shape_type = ellipse<2,T>;
volume (shape_type const &target):
m_target (target)
{ ; }
/// Generate a random point within the ellipse.
template <typename GeneratorT>
static auto
eval (ellipse<2,T> shape, GeneratorT &&g)
auto
eval (GeneratorT &&g) const noexcept
{
// We use a two step process:
// * Generate a point within a unit sphere
@ -177,9 +186,13 @@ namespace cruft::geom {
std::sin (phi)
};
auto const offset = circle_pos * rho * shape.radius;
return shape.origin + offset.template as<cruft::vector> ();
auto const offset = circle_pos * rho * m_target.radius;
return m_target.origin + offset.template as<cruft::vector> ();
}
private:
shape_type const &m_target;
};

2
geom/sample.cpp → geom/ops.cpp
View File

@ -6,4 +6,4 @@
* Copyright 2018 Danny Robson <danny@nerdcruft.net>
*/
#include "sample.hpp"
#include "ops.hpp"

18
geom/ops.hpp
View File

@ -6,8 +6,7 @@
* Copyright 2015-2018 Danny Robson <danny@nerdcruft.net>
*/
#ifndef __UTIL_GEOM_OPS_HPP
#define __UTIL_GEOM_OPS_HPP
#pragma once
#include "./fwd.hpp"
#include "../point.hpp"
@ -17,6 +16,9 @@
///////////////////////////////////////////////////////////////////////////////
namespace cruft::geom {
/// Tests whether any part of two shapes overlap.
///
/// Returns true even if the the edges of the shapes are co-linear.
template <
size_t S,
typename T,
@ -26,6 +28,8 @@ namespace cruft::geom {
bool
intersects (A<S,T>, B<S,T>);
/// Returns a minimum squared distance between two shapes.
template <
size_t S,
typename T,
@ -50,6 +54,8 @@ namespace cruft::geom {
T
distance (A<S,T>, B<S,T>);
/// Returns an AABB for the supplied shape.
template <
size_t S,
typename T,
@ -58,6 +64,8 @@ namespace cruft::geom {
aabb<S,T>
bounds (K<S,T>);
/// Returns a maximum distance across a shape.
template <
size_t S,
typename T,
@ -66,6 +74,8 @@ namespace cruft::geom {
T
diameter (K<S,T>);
/// Returns the closest point on a shape to the supplied point.
template <
size_t S,
typename T,
@ -74,6 +84,7 @@ namespace cruft::geom {
point<S,T>
closest (K<S,T>, point<S,T>);
template <
size_t S,
typename T,
@ -82,6 +93,7 @@ namespace cruft::geom {
vector<S,T>
magnitude (K<S,T>);
template <
size_t S,
typename T,
@ -90,5 +102,3 @@ namespace cruft::geom {
K<S,T>
scale (K<S,T>, T);
}
#endif

175
geom/sample.hpp
View File

@ -1,175 +0,0 @@
/*
* This Source Code Form is subject to the terms of the Mozilla Public
* License, v. 2.0. If a copy of the MPL was not distributed with this
* file, You can obtain one at http://mozilla.org/MPL/2.0/.
*
* Copyright 2015-2018 Danny Robson <danny@nerdcruft.net>
*/
#pragma once
#include "../coord/fwd.hpp"
#include "../extent.hpp"
#include "../random.hpp"
#include "ops.hpp"
#include <cstddef>
#include <random>
///////////////////////////////////////////////////////////////////////////////
namespace cruft::geom {
/// A function object that selects a uniformly random point inside a shape
/// using a provided random generator. The point will lie within the shape,
/// inclusive of boundaries.
///
/// May be specialised for arbitrary shapes but uses rejection sampling
/// as a safe default. This implies that execution may not take a constant
/// time.
template <typename ShapeT>
struct sampler {
/// Returns a point which lies within the supplied shape, inclusive
/// of borders.
template <typename GeneratorT>
static auto
eval (ShapeT const &shape, GeneratorT &&g)
{
auto b = bounds (shape);
while (true) {
auto p = sample (b, g);
if (intersects (shape, p))
return p;
}
}
};
///////////////////////////////////////////////////////////////////////////
/// A convenience function that calls sample::fn to select a random point
/// in a provided shape.
template <
typename ShapeT,
typename GeneratorT // random generator
>
auto
sample (ShapeT const &shape, GeneratorT &&gen)
{
return sampler<ShapeT>::eval (shape, std::forward<GeneratorT> (gen));
}
///////////////////////////////////////////////////////////////////////////
std::vector<cruft::point2f>
poisson_sample (cruft::extent2i, float distance, int samples);
///////////////////////////////////////////////////////////////////////////
namespace surface {
/// A generator of samples that lie on the surface of a shape
template <typename ShapeT>
class sampler;
template <typename ShapeT>
sampler (ShapeT const&) -> sampler<std::decay_t<ShapeT>>;
//---------------------------------------------------------------------
template <size_t S, typename T>
class sampler<cruft::extent<S,T>> {
public:
sampler (cruft::extent<S,T> _target):
target (_target)
{ ; }
template <typename GeneratorT>
cruft::point<S,T>
operator() (GeneratorT &&gen) const noexcept
{
cruft::point<S,T> p;
for (size_t i = 0; i < S; ++i)
p[i] = random::uniform (T{0}, target[i], gen);
return p;
}
private:
cruft::extent<S,T> target;
};
/// Approximate a poisson disc sampling through the "Mitchell's Best
/// Candidate" algorithm.
///
/// Try to keep adding a new point to a set. Each new point is the
/// best of a set of candidates. The 'best' is the point that is
/// furthest from all selected points.
///
/// \return A vector of the computed points
template <typename SamplerT, typename DistanceT, typename GeneratorT>
auto
poisson (SamplerT const &target,
GeneratorT &&gen,
DistanceT &&minimum_distance,
size_t candidates_count)
{
using point_type = decltype (target (gen));
using value_type = typename point_type::value_type;
std::vector<point_type> selected;
std::vector<point_type> candidates;
// prime the found elements list with an initial point we can
// perform distance calculations on
selected.push_back (target (gen));
// keep trying to add one more new point
while (1) {
// generate a group of candidate points
candidates.clear ();
std::generate_n (
std::back_inserter (candidates),
candidates_count,
[&] (void) {
return target (gen);
}
);
// find the point whose minimum distance to the existing
// points is the greatest (while also being greater than the
// required minimum distance);
auto best_distance2 = std::numeric_limits<value_type>::lowest ();
size_t best_index = 0;
for (size_t i = 0; i < candidates.size (); ++i) {
auto const p = candidates[i];
auto d2 = std::numeric_limits<value_type>::max ();
// find the minimum distance from this candidate to the
// selected points
for (auto q: selected)
d2 = cruft::min (d2, cruft::distance2 (p, q));
// record if it's the furthest away
if (d2 > best_distance2 && d2 > cruft::pow (minimum_distance (p), 2)) {
best_distance2 = d2;
best_index = i;
}
}
// if we didn't find a suitable point then we give up and
// return the points we found, otherwise record the best point
if (best_distance2 <= 0)
break;
selected.push_back (candidates[best_index]);
}
return selected;
}
}
}

33
geom/sample/fwd.hpp Normal file
View File

@ -0,0 +1,33 @@
/*
* This Source Code Form is subject to the terms of the Mozilla Public
* License, v. 2.0. If a copy of the MPL was not distributed with this
* file, You can obtain one at http://mozilla.org/MPL/2.0/.
*
* Copyright 2015-2018 Danny Robson <danny@nerdcruft.net>
*/
#pragma once
namespace cruft::geom::sample {
/// Provides an interface to randomly sample points within the volume of
/// a supplied shape.
///
/// This class _should_ be specialised for performance. Though a basic
/// rejection sampling implementation is defined as a fallback.
template <typename ShapeT>
class volume;
template <typename ShapeT>
volume (ShapeT const&) -> volume<std::decay_t<ShapeT>>;
/// Provides an interface to randomly sample points across the surface of
/// a supplied shape.
///
/// The class _must_ be specialised for each shape.
template <typename ShapeT>
class surface;
template <typename ShapeT>
surface (ShapeT const&) -> surface<std::decay_t<ShapeT>>;
}

105
geom/sample/surface.hpp Normal file
View File

@ -0,0 +1,105 @@
/*
* This Source Code Form is subject to the terms of the Mozilla Public
* License, v. 2.0. If a copy of the MPL was not distributed with this
* file, You can obtain one at http://mozilla.org/MPL/2.0/.
*
* Copyright 2015-2018 Danny Robson <danny@nerdcruft.net>
*/
#pragma once
#include "volume.hpp"
#include "../ops.hpp"
#include "../../coord/fwd.hpp"
#include "../../extent.hpp"
#include "../../random.hpp"
#include <cstddef>
#include <random>
///////////////////////////////////////////////////////////////////////////////
namespace cruft::geom::sample {
/// Approximate a poisson disc sampling through the "Mitchell's Best
/// Candidate" algorithm.
///
/// Try to keep adding a new point to a set. Each new point is the
/// best of a set of candidates. The 'best' is the point that is
/// furthest from all selected points.
///
/// \return A vector of the computed points
template <typename SamplerT, typename DistanceT, typename GeneratorT>
auto
poisson (SamplerT const &target,
GeneratorT &&gen,
DistanceT &&minimum_distance,
size_t candidates_count)
{
using point_type = decltype (target.eval (gen));
using value_type = typename point_type::value_type;
std::vector<point_type> selected;
std::vector<point_type> candidates;
// prime the found elements list with an initial point we can
// perform distance calculations on
selected.push_back (target.eval (gen));
// keep trying to add one more new point
while (1) {
// generate a group of candidate points
candidates.clear ();
std::generate_n (
std::back_inserter (candidates),
candidates_count,
[&] (void) {
return target.eval (gen);
}
);
// find the point whose minimum distance to the existing
// points is the greatest (while also being greater than the
// required minimum distance);
auto best_distance2 = std::numeric_limits<value_type>::lowest ();
size_t best_index = 0;
for (size_t i = 0; i < candidates.size (); ++i) {
auto const p = candidates[i];
auto d2 = std::numeric_limits<value_type>::max ();
// find the minimum distance from this candidate to the
// selected points
for (auto q: selected)
d2 = cruft::min (d2, cruft::distance2 (p, q));
// record if it's the furthest away
if (d2 > best_distance2 && d2 > cruft::pow (minimum_distance (p), 2)) {
best_distance2 = d2;
best_index = i;
}
}
// if we didn't find a suitable point then we give up and
// return the points we found, otherwise record the best point
if (best_distance2 <= 0)
break;
selected.push_back (candidates[best_index]);
}
return selected;
}
/// A surface sampler specialisation for 2d extents.
///
/// The actual work is handed off to the volume sampler, as it's
/// equivalent in 2 dimensions.
template <typename T>
class surface<extent<2,T>> :
public volume<extent<2,T>>
{ using volume<extent<2,T>>::volume; };
}

105
geom/sample/volume.hpp Normal file
View File

@ -0,0 +1,105 @@
/*
* This Source Code Form is subject to the terms of the Mozilla Public
* License, v. 2.0. If a copy of the MPL was not distributed with this
* file, You can obtain one at http://mozilla.org/MPL/2.0/.
*
* Copyright 2015-2018 Danny Robson <danny@nerdcruft.net>
*/
#pragma once
#include "../../coord/fwd.hpp"
#include "../../extent.hpp"
#include "../../random.hpp"
#include "../ops.hpp"
#include <cstddef>
#include <random>
///////////////////////////////////////////////////////////////////////////////
namespace cruft::geom::sample {
///////////////////////////////////////////////////////////////////////
/// A convenience function that calls sample::fn to select a random
/// point in a provided shape.
///
/// This function is useful because we can perform use partial
/// specialisation if we pass the functionality off to a well known
/// class.
template <
typename ShapeT,
typename GeneratorT // random generator
>
auto
eval (ShapeT const &shape, GeneratorT &&gen)
{
volume query (shape);
return query.eval (std::forward<GeneratorT> (gen));
}
/// A function object that selects a uniformly random point inside a
/// shape using a provided random generator. The point will lie
/// within the shape, inclusive of boundaries.
///
/// This should be specialised for anything complex, but will work
/// (with a varying performance hit) for anything that provides
/// bounds` and `intersects`.
template <typename ShapeT>
class volume {
public:
explicit volume (ShapeT const &target):
m_target (target)
{ ; }
/// Returns a point which lies within the supplied shape, inclusive
/// of borders.
///
template <typename GeneratorT>
auto
eval (GeneratorT &&g) const
{
auto b = bounds (m_target);
while (true) {
auto p = ::cruft::geom::sample::eval (b, g);
if (intersects (m_target, p))
return p;
}
}
private:
ShapeT const &m_target;
};
///////////////////////////////////////////////////////////////////////
/// A specialisation of the volume sampler for extents.
template <size_t S, typename T>
class volume<cruft::extent<S,T>> {
public:
using shape_type = extent<S,T>;
explicit volume (shape_type const &target):
m_target (target)
{ ; }
template <typename GeneratorT>
auto
eval (GeneratorT &&gen) const noexcept
{
cruft::point<S,T> p;
for (size_t i = 0; i < S; ++i)
p[i] = random::uniform (T{0}, m_target[i], gen);
return p;
}
private:
shape_type const &m_target;
};
}

13
geom/tri.hpp
View File

@ -10,7 +10,7 @@
#define CRUFT_GEOM_TRI_HPP
#include "../point.hpp"
#include "sample.hpp"
#include "sample/fwd.hpp"
#include <cstddef>
#include <random>
@ -33,11 +33,13 @@ namespace cruft::geom {
using tri3f = tri<3,float>;
namespace surface {
namespace sample {
template <size_t S, typename T>
class sampler<tri<S,T>> {
class surface<tri<S,T>> {
public:
sampler (tri<S,T> _target):
using shape_type = tri<S,T>;
explicit surface (tri<S,T> _target):
base (_target.a),
v0 (_target.b - _target.a),
v1 (_target.c - _target.a)
@ -53,7 +55,8 @@ namespace cruft::geom {
private:
cruft::point<S,T> base;
cruft::vector<S,T> v0, v1;
cruft::vector<S,T> v0;
cruft::vector<S,T> v1;
};
};

11
tools/poisson.cpp
View File

@ -1,7 +1,7 @@
#include "cmdopt.hpp"
#include "functor.hpp"
#include "geom/aabb.hpp"
#include "geom/sample.hpp"
#include "geom/sample/surface.hpp"
#include <cstdlib>
#include <iostream>
@ -25,10 +25,11 @@ main (int argc, char **argv)
opts.scan (argc, argv);
std::cout << "<svg height='" << area.h << "' width='" << area.h << "'>";
for (auto p: cruft::geom::surface::poisson (cruft::geom::surface::sampler (area),
std::default_random_engine (),
cruft::functor::constant (distance),
samples))
namespace sample = cruft::geom::sample;
for (auto const &p: sample::poisson (sample::surface {area},
std::default_random_engine (),
cruft::functor::constant (distance),
samples))
{
std::cout << "<circle cx='" << p.x << "' cy='" << p.y << "' r='1' />";
}