Skip to content

Commit 361861c

Browse files
tinko92tinko92
committed
Add some specialized distributions for rings that are triangles (cartesian 2d only at this time) or convex and a general rejection sampling strategy based on convex hull.
1 parent e981c46 commit 361861c

File tree

3 files changed

+259
-0
lines changed

3 files changed

+259
-0
lines changed

‎include/boost/geometry/extensions/random/strategies/agnostic/uniform_convex_fan.hpp‎

Lines changed: 95 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,95 @@
1+
// Boost.Geometry (aka GGL, Generic Geometry Library)
2+
3+
// Copyright (c) 2019 Tinko Bartels, Berlin, Germany.
4+
5+
// Use, modification and distribution is subject to the Boost Software License,
6+
// Version 1.0. (See accompanying file LICENSE_1_0.txt or copy at
7+
// http://www.boost.org/LICENSE_1_0.txt)
8+
9+
#ifndef BOOST_GEOMETRY_EXTENSIONS_RANDOM_STRATEGIES_AGNOSTIC_UNIFORM_CONVEX_FAN_HPP
10+
#define BOOST_GEOMETRY_EXTENSIONS_RANDOM_STRATEGIES_AGNOSTIC_UNIFORM_CONVEX_FAN_HPP
11+
12+
#include <random>
13+
#include <vector>
14+
#include <cstdlib>
15+
#include <cmath>
16+
17+
#include <boost/geometry/algorithms/equals.hpp>
18+
#include <boost/geometry/algorithms/transform.hpp>
19+
#include <boost/geometry/algorithms/within.hpp>
20+
21+
namespace boost { namespace geometry
22+
{
23+
24+
namespace strategy { namespace uniform_point_distribution {
25+
26+
// The following strategy is suitable for convex rings and polygons with
27+
// non-empty interior.
28+
template
29+
<
30+
typename Point,
31+
typename DomainGeometry,
32+
typename TriangleStrategy,
33+
typename SideStrategy //Actually, we need a triangle area strategy here.
34+
>
35+
struct uniform_convex_fan
36+
{
37+
private:
38+
std::vector<double> accumulated_areas;
39+
// It is hard to see a reason not to use double here. If a triangles
40+
// relative size is smaller than doubles epsilon, it is too unlikely to
41+
// realistically occur in a random sample anyway.
42+
public:
43+
uniform_convex_fan(DomainGeometry const& g)
44+
{
45+
accumulated_areas.push_back(0);
46+
for (int i = 2 ; i < g.size() ; ++i) {
47+
accumulated_areas.push_back(
48+
accumulated_areas.back() +
49+
std::abs(SideStrategy::template side_value<double, double>(
50+
*g.begin(),
51+
*(g.begin() + i - 1),
52+
*(g.begin() + i))));
53+
}
54+
}
55+
bool equals(DomainGeometry const& l_domain,
56+
DomainGeometry const& r_domain,
57+
uniform_convex_fan const& r_strategy) const
58+
{
59+
if( l_domain.size() != r_domain.size() ) return false;
60+
for (int i = 0; i < l_domain.size(); ++i) {
61+
if( !boost::geometry::equals(*(l_domain.begin() + i),
62+
*(r_domain.begin() + i)))
63+
return false;
64+
}
65+
return true;
66+
}
67+
template<typename Gen>
68+
Point apply(Gen& g, DomainGeometry const& d)
69+
{
70+
std::uniform_real_distribution<double> dist(0, 1);
71+
double r = dist(g) * accumulated_areas.back(),
72+
s = dist(g);
73+
std::size_t i = std::distance(
74+
accumulated_areas.begin(),
75+
std::lower_bound(accumulated_areas.begin(),
76+
accumulated_areas.end(),
77+
r));
78+
return TriangleStrategy::template map
79+
<
80+
double
81+
>(* d.begin(),
82+
*(d.begin() + i),
83+
*(d.begin() + i + 1),
84+
( r - accumulated_areas[ i - 1 ]) /
85+
( accumulated_areas[ i ] - accumulated_areas[ i - 1 ] ),
86+
s);
87+
}
88+
void reset(DomainGeometry const&) {};
89+
};
90+
91+
}} // namespace strategy::uniform_point_distribution
92+
93+
}} // namespace boost::geometry
94+
95+
#endif // BOOST_GEOMETRY_EXTENSIONS_RANDOM_STRATEGIES_AGNOSTIC_UNIFORM_CONVEX_FAN_HPP

‎include/boost/geometry/extensions/random/strategies/agnostic/uniform_convex_hull_rejection.hpp‎

Lines changed: 82 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,82 @@
1+
// Boost.Geometry (aka GGL, Generic Geometry Library)
2+
3+
// Copyright (c) 2019 Tinko Bartels, Berlin, Germany.
4+
5+
// Use, modification and distribution is subject to the Boost Software License,
6+
// Version 1.0. (See accompanying file LICENSE_1_0.txt or copy at
7+
// http://www.boost.org/LICENSE_1_0.txt)
8+
9+
#ifndef BOOST_GEOMETRY_EXTENSIONS_RANDOM_STRATEGIES_AGNOSTIC_UNIFORM_CONVEX_HULL_REJECTION_HPP
10+
#define BOOST_GEOMETRY_EXTENSIONS_RANDOM_STRATEGIES_AGNOSTIC_UNIFORM_CONVEX_HULL_REJECTION_HPP
11+
12+
#include <random>
13+
#include <vector>
14+
15+
#include <boost/geometry/algorithms/equals.hpp>
16+
#include <boost/geometry/algorithms/transform.hpp>
17+
#include <boost/geometry/algorithms/within.hpp>
18+
#include <boost/geometry/geometries/polygon.hpp>
19+
#include <boost/geometry/algorithms/convex_hull.hpp>
20+
21+
#include <boost/geometry/extensions/random/strategies/agnostic/uniform_convex_fan.hpp>
22+
23+
namespace boost { namespace geometry
24+
{
25+
26+
namespace strategy { namespace uniform_point_distribution {
27+
28+
// The following strategy is suitable for geometries for which
29+
// a Triangle sampling strategy can be provided
30+
// a SideStrategy that compues the triangle area can be provided.
31+
template
32+
<
33+
typename Point,
34+
typename DomainGeometry,
35+
typename TriangleStrategy,
36+
typename SideStrategy //Actually, we need a triangle area strategy here.
37+
>
38+
struct uniform_convex_hull_rejection
39+
{
40+
private:
41+
typedef typename point_type<DomainGeometry>::type domain_point_type;
42+
typedef boost::geometry::model::ring<domain_point_type> ring;
43+
ring hull;
44+
uniform_convex_fan<Point, ring, TriangleStrategy, SideStrategy>
45+
m_strategy;
46+
public:
47+
uniform_convex_hull_rejection(DomainGeometry const& g) : m_strategy(hull)
48+
{
49+
boost::geometry::convex_hull(g, hull);
50+
m_strategy =
51+
uniform_convex_fan
52+
<
53+
Point,
54+
ring,
55+
TriangleStrategy,
56+
SideStrategy
57+
>(hull);
58+
}
59+
bool equals(DomainGeometry const& l_domain,
60+
DomainGeometry const& r_domain,
61+
uniform_convex_hull_rejection const& r_strategy) const
62+
{
63+
return boost::geometry::equals(l_domain, r_domain)
64+
&& m_strategy.equals(l_domain, r_domain, r_strategy.m_strategy);
65+
}
66+
template<typename Gen>
67+
Point apply(Gen& g, DomainGeometry const& d)
68+
{
69+
Point p;
70+
do{
71+
p = m_strategy.apply(g, hull);
72+
}while( !boost::geometry::within(p, d) );
73+
return p;
74+
}
75+
void reset(DomainGeometry const&) {};
76+
};
77+
78+
}} // namespace strategy::uniform_point_distribution
79+
80+
}} // namespace boost::geometry
81+
82+
#endif // BOOST_GEOMETRY_EXTENSIONS_RANDOM_STRATEGIES_AGNOSTIC_UNIFORM_CONVEX_HULL_REJECTION_HPP

‎include/boost/geometry/extensions/random/strategies/cartesian/uniform_point_distribution_triangle.hpp‎

Lines changed: 82 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,82 @@
1+
// Boost.Geometry (aka GGL, Generic Geometry Library)
2+
3+
// Copyright (c) 2019 Tinko Bartels, Berlin, Germany.
4+
5+
// Use, modification and distribution is subject to the Boost Software License,
6+
// Version 1.0. (See accompanying file LICENSE_1_0.txt or copy at
7+
// http://www.boost.org/LICENSE_1_0.txt)
8+
9+
#ifndef BOOST_GEOMETRY_EXTENSIONS_RANDOM_STRATEGIES_CARTESIAN_UNIFORM_POINT_DISTRIBUTION_TRIANGLE_HPP
10+
#define BOOST_GEOMETRY_EXTENSIONS_RANDOM_STRATEGIES_CARTESIAN_UNIFORM_POINT_DISTRIBUTION_TRIANGLE_HPP
11+
12+
#include <random>
13+
#include <cmath>
14+
15+
#include <boost/geometry/arithmetic/arithmetic.hpp>
16+
#include <boost/geometry/algorithms/assign.hpp>
17+
#include <boost/geometry/core/access.hpp>
18+
#include <boost/geometry/core/point_type.hpp>
19+
20+
namespace boost { namespace geometry
21+
{
22+
23+
namespace strategy { namespace uniform_point_distribution {
24+
25+
//The following strategy is suitable for rings or polygons of three points.
26+
template
27+
<
28+
typename Point,
29+
typename DomainGeometry
30+
>
31+
struct uniform_2d_cartesian_triangle
32+
{
33+
private:
34+
typedef typename point_type<DomainGeometry>::type domain_point_type;
35+
public:
36+
uniform_2d_cartesian_triangle(DomainGeometry const& g) {}
37+
bool equals(DomainGeometry const& l_domain,
38+
DomainGeometry const& r_domain,
39+
uniform_2d_cartesian_triangle const& r_strategy) const
40+
{
41+
return boost::geometry::equals(l_domain.domain(), r_domain.domain());
42+
}
43+
44+
template<typename sample_type>
45+
static Point map(domain_point_type const& p1,
46+
domain_point_type const& p2,
47+
domain_point_type const& p3,
48+
sample_type const& s1,
49+
sample_type const& s2)
50+
{
51+
Point out;
52+
sample_type r1 = std::sqrt(s1);
53+
set<0>(out, (1 - r1) * get<0>(p1)
54+
+ ( r1 * (1 - s2) ) * get<0>(p2)
55+
+ ( s2 * r1 * get<0>(p3)));
56+
set<1>(out, (1 - r1) * get<1>(p1)
57+
+ ( r1 * (1 - s2) ) * get<1>(p2)
58+
+ ( s2 * r1 * get<1>(p3)));
59+
return out;
60+
}
61+
62+
template<typename Gen>
63+
Point apply(Gen& g, DomainGeometry const& d)
64+
{
65+
typedef typename select_most_precise
66+
<
67+
typename coordinate_type<DomainGeometry>::type,
68+
double
69+
>::type sample_type;
70+
std::uniform_real_distribution<sample_type> real_dist(0, 1);
71+
sample_type s1 = real_dist(g),
72+
s2 = real_dist(g);
73+
return map(*d.begin(), *(d.begin() + 1), *(d.begin() + 2), s1, s2);
74+
}
75+
void reset(DomainGeometry const&) {};
76+
};
77+
78+
}} // namespace strategy::uniform_point_distribution
79+
80+
}} // namespace boost::geometry
81+
82+
#endif // BOOST_GEOMETRY_EXTENSIONS_RANDOM_STRATEGIES_CARTESIAN_UNIFORM_POINT_DISTRIBUTION_TRIANGLE_HPP

0 commit comments

Comments
 (0)