12 #ifndef RANDOM_POINT_GENERATORS_H_
13 #define RANDOM_POINT_GENERATORS_H_
15 #include <CGAL/number_utils.h>
16 #include <CGAL/Random.h>
17 #include <CGAL/point_generators_d.h>
18 #include <CGAL/version.h>
23 #if CGAL_VERSION_NR < 1041101000
24 # error random_point_generators is only available for CGAL >= 4.11
35 template <
typename Kernel>
36 typename Kernel::Point_d construct_point(
const Kernel &k,
37 typename Kernel::FT x1,
typename Kernel::FT x2) {
38 typename Kernel::FT tab[2];
41 return k.construct_point_d_object()(2, &tab[0], &tab[2]);
46 template <
typename Kernel>
47 typename Kernel::Point_d construct_point(
const Kernel &k,
48 typename Kernel::FT x1,
typename Kernel::FT x2,
typename Kernel::FT x3) {
49 typename Kernel::FT tab[3];
53 return k.construct_point_d_object()(3, &tab[0], &tab[3]);
58 template <
typename Kernel>
59 typename Kernel::Point_d construct_point(
const Kernel &k,
60 typename Kernel::FT x1,
typename Kernel::FT x2,
typename Kernel::FT x3,
61 typename Kernel::FT x4) {
62 typename Kernel::FT tab[4];
67 return k.construct_point_d_object()(4, &tab[0], &tab[4]);
72 template <
typename Kernel>
73 typename Kernel::Point_d construct_point(
const Kernel &k,
74 typename Kernel::FT x1,
typename Kernel::FT x2,
typename Kernel::FT x3,
75 typename Kernel::FT x4,
typename Kernel::FT x5) {
76 typename Kernel::FT tab[5];
82 return k.construct_point_d_object()(5, &tab[0], &tab[5]);
87 template <
typename Kernel>
88 typename Kernel::Point_d construct_point(
const Kernel &k,
89 typename Kernel::FT x1,
typename Kernel::FT x2,
typename Kernel::FT x3,
90 typename Kernel::FT x4,
typename Kernel::FT x5,
typename Kernel::FT x6) {
91 typename Kernel::FT tab[6];
98 return k.construct_point_d_object()(6, &tab[0], &tab[6]);
101 template <
typename Kernel>
102 std::vector<typename Kernel::Point_d> generate_points_on_plane(std::size_t num_points,
int intrinsic_dim,
104 double coord_min = -5.,
double coord_max = 5.) {
105 typedef typename Kernel::Point_d Point;
106 typedef typename Kernel::FT FT;
109 std::vector<Point> points;
110 points.reserve(num_points);
111 for (std::size_t i = 0; i < num_points;) {
112 std::vector<FT> pt(ambient_dim, FT(0));
113 for (
int j = 0; j < intrinsic_dim; ++j)
114 pt[j] = rng.get_double(coord_min, coord_max);
116 Point p = k.construct_point_d_object()(ambient_dim, pt.begin(), pt.end());
123 template <
typename Kernel>
124 std::vector<typename Kernel::Point_d> generate_points_on_moment_curve(std::size_t num_points,
int dim,
125 typename Kernel::FT min_x,
126 typename Kernel::FT max_x) {
127 typedef typename Kernel::Point_d Point;
128 typedef typename Kernel::FT FT;
131 std::vector<Point> points;
132 points.reserve(num_points);
133 for (std::size_t i = 0; i < num_points;) {
134 FT x = rng.get_double(min_x, max_x);
135 std::vector<FT> coords;
137 for (
int p = 1; p <= dim; ++p)
138 coords.push_back(std::pow(CGAL::to_double(x), p));
139 Point p = k.construct_point_d_object()(
140 dim, coords.begin(), coords.end());
149 template <
typename Kernel>
150 std::vector<typename Kernel::Point_d> generate_points_on_torus_3D(std::size_t num_points,
double R,
double r,
151 bool uniform =
false) {
152 typedef typename Kernel::Point_d Point;
153 typedef typename Kernel::FT FT;
158 std::size_t num_lines = (std::size_t)sqrt(num_points);
160 std::vector<Point> points;
161 points.reserve(num_points);
162 for (std::size_t i = 0; i < num_points;) {
165 std::size_t k1 = i / num_lines;
166 std::size_t k2 = i % num_lines;
167 u = 6.2832 * k1 / num_lines;
168 v = 6.2832 * k2 / num_lines;
170 u = rng.get_double(0, 6.2832);
171 v = rng.get_double(0, 6.2832);
173 Point p = construct_point(k,
174 (R + r * std::cos(u)) * std::cos(v),
175 (R + r * std::cos(u)) * std::sin(v),
184 template <
typename Kernel,
typename OutputIterator>
185 static void generate_uniform_points_on_torus_d(
const Kernel &k,
int dim, std::size_t num_slices,
187 double radius_noise_percentage = 0.,
188 std::vector<typename Kernel::FT> current_point =
189 std::vector<typename Kernel::FT>()) {
191 int point_size =
static_cast<int>(current_point.size());
192 if (point_size == 2 * dim) {
193 *out++ = k.construct_point_d_object()(point_size, current_point.begin(), current_point.end());
195 for (std::size_t slice_idx = 0; slice_idx < num_slices; ++slice_idx) {
196 double radius_noise_ratio = 1.;
197 if (radius_noise_percentage > 0.) {
198 radius_noise_ratio = rng.get_double(
199 (100. - radius_noise_percentage) / 100.,
200 (100. + radius_noise_percentage) / 100.);
202 std::vector<typename Kernel::FT> cp2 = current_point;
203 double alpha = 6.2832 * slice_idx / num_slices;
204 cp2.push_back(radius_noise_ratio * std::cos(alpha));
205 cp2.push_back(radius_noise_ratio * std::sin(alpha));
206 generate_uniform_points_on_torus_d(
207 k, dim, num_slices, out, radius_noise_percentage, cp2);
212 template <
typename Kernel>
213 std::vector<typename Kernel::Point_d> generate_points_on_torus_d(std::size_t num_points,
int dim,
bool uniform =
false,
214 double radius_noise_percentage = 0.) {
215 typedef typename Kernel::Point_d Point;
216 typedef typename Kernel::FT FT;
220 std::vector<Point> points;
221 points.reserve(num_points);
223 std::size_t num_slices = (std::size_t)std::pow(num_points, 1. / dim);
224 generate_uniform_points_on_torus_d(
225 k, dim, num_slices, std::back_inserter(points), radius_noise_percentage);
227 for (std::size_t i = 0; i < num_points;) {
228 double radius_noise_ratio = 1.;
229 if (radius_noise_percentage > 0.) {
230 radius_noise_ratio = rng.get_double(
231 (100. - radius_noise_percentage) / 100.,
232 (100. + radius_noise_percentage) / 100.);
234 std::vector<typename Kernel::FT> pt;
236 for (
int curdim = 0; curdim < dim; ++curdim) {
237 FT alpha = rng.get_double(0, 6.2832);
238 pt.push_back(radius_noise_ratio * std::cos(alpha));
239 pt.push_back(radius_noise_ratio * std::sin(alpha));
242 Point p = k.construct_point_d_object()(pt.begin(), pt.end());
250 template <
typename Kernel>
251 std::vector<typename Kernel::Point_d> generate_points_on_sphere_d(std::size_t num_points,
int dim,
double radius,
252 double radius_noise_percentage = 0.) {
253 typedef typename Kernel::Point_d Point;
256 CGAL::Random_points_on_sphere_d<Point> generator(dim, radius);
257 std::vector<Point> points;
258 points.reserve(num_points);
259 for (std::size_t i = 0; i < num_points;) {
260 Point p = *generator++;
261 if (radius_noise_percentage > 0.) {
262 double radius_noise_ratio = rng.get_double(
263 (100. - radius_noise_percentage) / 100.,
264 (100. + radius_noise_percentage) / 100.);
266 typename Kernel::Point_to_vector_d k_pt_to_vec =
267 k.point_to_vector_d_object();
268 typename Kernel::Vector_to_point_d k_vec_to_pt =
269 k.vector_to_point_d_object();
270 typename Kernel::Scaled_vector_d k_scaled_vec =
271 k.scaled_vector_d_object();
272 p = k_vec_to_pt(k_scaled_vec(k_pt_to_vec(p), radius_noise_ratio));
280 template <
typename Kernel>
281 std::vector<typename Kernel::Point_d> generate_points_in_ball_d(std::size_t num_points,
int dim,
double radius) {
282 typedef typename Kernel::Point_d Point;
285 CGAL::Random_points_in_ball_d<Point> generator(dim, radius);
286 std::vector<Point> points;
287 points.reserve(num_points);
288 for (std::size_t i = 0; i < num_points;) {
289 Point p = *generator++;
296 template <
typename Kernel>
297 std::vector<typename Kernel::Point_d> generate_points_in_cube_d(std::size_t num_points,
int dim,
double radius) {
298 typedef typename Kernel::Point_d Point;
301 CGAL::Random_points_in_cube_d<Point> generator(dim, radius);
302 std::vector<Point> points;
303 points.reserve(num_points);
304 for (std::size_t i = 0; i < num_points;) {
305 Point p = *generator++;
312 template <
typename Kernel>
313 std::vector<typename Kernel::Point_d> generate_points_on_two_spheres_d(std::size_t num_points,
int dim,
double radius,
314 double distance_between_centers,
315 double radius_noise_percentage = 0.) {
316 typedef typename Kernel::FT FT;
317 typedef typename Kernel::Point_d Point;
318 typedef typename Kernel::Vector_d Vector;
321 CGAL::Random_points_on_sphere_d<Point> generator(dim, radius);
322 std::vector<Point> points;
323 points.reserve(num_points);
325 std::vector<FT> t(dim, FT(0));
326 t[0] = distance_between_centers;
327 Vector c1_to_c2(t.begin(), t.end());
329 for (std::size_t i = 0; i < num_points;) {
330 Point p = *generator++;
331 if (radius_noise_percentage > 0.) {
332 double radius_noise_ratio = rng.get_double(
333 (100. - radius_noise_percentage) / 100.,
334 (100. + radius_noise_percentage) / 100.);
336 typename Kernel::Point_to_vector_d k_pt_to_vec =
337 k.point_to_vector_d_object();
338 typename Kernel::Vector_to_point_d k_vec_to_pt =
339 k.vector_to_point_d_object();
340 typename Kernel::Scaled_vector_d k_scaled_vec =
341 k.scaled_vector_d_object();
342 p = k_vec_to_pt(k_scaled_vec(k_pt_to_vec(p), radius_noise_ratio));
345 typename Kernel::Translated_point_d k_transl =
346 k.translated_point_d_object();
347 Point p2 = k_transl(p, c1_to_c2);
349 points.push_back(p2);
357 template <
typename Kernel>
358 std::vector<typename Kernel::Point_d> generate_points_on_3sphere_and_circle(std::size_t num_points,
359 double sphere_radius) {
360 typedef typename Kernel::FT FT;
361 typedef typename Kernel::Point_d Point;
364 CGAL::Random_points_on_sphere_d<Point> generator(3, sphere_radius);
365 std::vector<Point> points;
366 points.reserve(num_points);
368 typename Kernel::Compute_coordinate_d k_coord =
369 k.compute_coordinate_d_object();
370 for (std::size_t i = 0; i < num_points;) {
371 Point p_sphere = *generator++;
373 FT alpha = rng.get_double(0, 6.2832);
374 std::vector<FT> pt(5);
375 pt[0] = k_coord(p_sphere, 0);
376 pt[1] = k_coord(p_sphere, 1);
377 pt[2] = k_coord(p_sphere, 2);
378 pt[3] = std::cos(alpha);
379 pt[4] = std::sin(alpha);
380 Point p(pt.begin(), pt.end());
388 template <
typename Kernel>
389 std::vector<typename Kernel::Point_d> generate_points_on_klein_bottle_3D(std::size_t num_points,
double a,
double b,
390 bool uniform =
false) {
391 typedef typename Kernel::Point_d Point;
392 typedef typename Kernel::FT FT;
397 std::size_t num_lines = (std::size_t)sqrt(num_points);
399 std::vector<Point> points;
400 points.reserve(num_points);
401 for (std::size_t i = 0; i < num_points;) {
404 std::size_t k1 = i / num_lines;
405 std::size_t k2 = i % num_lines;
406 u = 6.2832 * k1 / num_lines;
407 v = 6.2832 * k2 / num_lines;
409 u = rng.get_double(0, 6.2832);
410 v = rng.get_double(0, 6.2832);
412 double tmp = cos(u / 2) * sin(v) - sin(u / 2) * sin(2. * v);
413 Point p = construct_point(k,
414 (a + b * tmp) * cos(u),
415 (a + b * tmp) * sin(u),
416 b * (sin(u / 2) * sin(v) + cos(u / 2) * sin(2. * v)));
424 template <
typename Kernel>
425 std::vector<typename Kernel::Point_d> generate_points_on_klein_bottle_4D(std::size_t num_points,
double a,
double b,
426 double noise = 0.,
bool uniform =
false) {
427 typedef typename Kernel::Point_d Point;
428 typedef typename Kernel::FT FT;
433 std::size_t num_lines = (std::size_t)sqrt(num_points);
435 std::vector<Point> points;
436 points.reserve(num_points);
437 for (std::size_t i = 0; i < num_points;) {
440 std::size_t k1 = i / num_lines;
441 std::size_t k2 = i % num_lines;
442 u = 6.2832 * k1 / num_lines;
443 v = 6.2832 * k2 / num_lines;
445 u = rng.get_double(0, 6.2832);
446 v = rng.get_double(0, 6.2832);
448 Point p = construct_point(k,
449 (a + b * cos(v)) * cos(u) + (noise == 0. ? 0. : rng.get_double(0, noise)),
450 (a + b * cos(v)) * sin(u) + (noise == 0. ? 0. : rng.get_double(0, noise)),
451 b * sin(v) * cos(u / 2) + (noise == 0. ? 0. : rng.get_double(0, noise)),
452 b * sin(v) * sin(u / 2) + (noise == 0. ? 0. : rng.get_double(0, noise)));
462 template <
typename Kernel>
463 std::vector<typename Kernel::Point_d>
464 generate_points_on_klein_bottle_variant_5D(
465 std::size_t num_points,
double a,
double b,
bool uniform =
false) {
466 typedef typename Kernel::Point_d Point;
467 typedef typename Kernel::FT FT;
472 std::size_t num_lines = (std::size_t)sqrt(num_points);
474 std::vector<Point> points;
475 points.reserve(num_points);
476 for (std::size_t i = 0; i < num_points;) {
479 std::size_t k1 = i / num_lines;
480 std::size_t k2 = i % num_lines;
481 u = 6.2832 * k1 / num_lines;
482 v = 6.2832 * k2 / num_lines;
484 u = rng.get_double(0, 6.2832);
485 v = rng.get_double(0, 6.2832);
487 FT x1 = (a + b * cos(v)) * cos(u);
488 FT x2 = (a + b * cos(v)) * sin(u);
489 FT x3 = b * sin(v) * cos(u / 2);
490 FT x4 = b * sin(v) * sin(u / 2);
491 FT x5 = x1 + x2 + x3 + x4;
493 Point p = construct_point(k, x1, x2, x3, x4, x5);
502 #endif // RANDOM_POINT_GENERATORS_H_