23 #ifndef RANDOM_POINT_GENERATORS_H_ 24 #define RANDOM_POINT_GENERATORS_H_ 26 #include <CGAL/number_utils.h> 27 #include <CGAL/Random.h> 28 #include <CGAL/point_generators_d.h> 40 template <
typename Kernel>
41 typename Kernel::Point_d construct_point(
const Kernel &k,
42 typename Kernel::FT x1,
typename Kernel::FT x2) {
43 typename Kernel::FT tab[2];
46 return k.construct_point_d_object()(2, &tab[0], &tab[2]);
51 template <
typename Kernel>
52 typename Kernel::Point_d construct_point(
const Kernel &k,
53 typename Kernel::FT x1,
typename Kernel::FT x2,
typename Kernel::FT x3) {
54 typename Kernel::FT tab[3];
58 return k.construct_point_d_object()(3, &tab[0], &tab[3]);
63 template <
typename Kernel>
64 typename Kernel::Point_d construct_point(
const Kernel &k,
65 typename Kernel::FT x1,
typename Kernel::FT x2,
typename Kernel::FT x3,
66 typename Kernel::FT x4) {
67 typename Kernel::FT tab[4];
72 return k.construct_point_d_object()(4, &tab[0], &tab[4]);
77 template <
typename Kernel>
78 typename Kernel::Point_d construct_point(
const Kernel &k,
79 typename Kernel::FT x1,
typename Kernel::FT x2,
typename Kernel::FT x3,
80 typename Kernel::FT x4,
typename Kernel::FT x5) {
81 typename Kernel::FT tab[5];
87 return k.construct_point_d_object()(5, &tab[0], &tab[5]);
92 template <
typename Kernel>
93 typename Kernel::Point_d construct_point(
const Kernel &k,
94 typename Kernel::FT x1,
typename Kernel::FT x2,
typename Kernel::FT x3,
95 typename Kernel::FT x4,
typename Kernel::FT x5,
typename Kernel::FT x6) {
96 typename Kernel::FT tab[6];
103 return k.construct_point_d_object()(6, &tab[0], &tab[6]);
106 template <
typename Kernel>
107 std::vector<typename Kernel::Point_d> generate_points_on_plane(std::size_t num_points,
int intrinsic_dim,
109 double coord_min = -5.,
double coord_max = 5.) {
110 typedef typename Kernel::Point_d Point;
111 typedef typename Kernel::FT FT;
114 std::vector<Point> points;
115 points.reserve(num_points);
116 for (std::size_t i = 0; i < num_points;) {
117 std::vector<FT> pt(ambient_dim, FT(0));
118 for (
int j = 0; j < intrinsic_dim; ++j)
119 pt[j] = rng.get_double(coord_min, coord_max);
121 Point p = k.construct_point_d_object()(ambient_dim, pt.begin(), pt.end());
128 template <
typename Kernel>
129 std::vector<typename Kernel::Point_d> generate_points_on_moment_curve(std::size_t num_points,
int dim,
130 typename Kernel::FT min_x,
131 typename Kernel::FT max_x) {
132 typedef typename Kernel::Point_d Point;
133 typedef typename Kernel::FT FT;
136 std::vector<Point> points;
137 points.reserve(num_points);
138 for (std::size_t i = 0; i < num_points;) {
139 FT x = rng.get_double(min_x, max_x);
140 std::vector<FT> coords;
142 for (
int p = 1; p <= dim; ++p)
143 coords.push_back(std::pow(CGAL::to_double(x), p));
144 Point p = k.construct_point_d_object()(
145 dim, coords.begin(), coords.end());
154 template <
typename Kernel>
155 std::vector<typename Kernel::Point_d> generate_points_on_torus_3D(std::size_t num_points,
double R,
double r,
156 bool uniform =
false) {
157 typedef typename Kernel::Point_d Point;
158 typedef typename Kernel::FT FT;
163 std::size_t num_lines = (std::size_t)sqrt(num_points);
165 std::vector<Point> points;
166 points.reserve(num_points);
167 for (std::size_t i = 0; i < num_points;) {
170 std::size_t k1 = i / num_lines;
171 std::size_t k2 = i % num_lines;
172 u = 6.2832 * k1 / num_lines;
173 v = 6.2832 * k2 / num_lines;
175 u = rng.get_double(0, 6.2832);
176 v = rng.get_double(0, 6.2832);
178 Point p = construct_point(k,
179 (R + r * std::cos(u)) * std::cos(v),
180 (R + r * std::cos(u)) * std::sin(v),
189 template <
typename Kernel,
typename OutputIterator>
190 static void generate_uniform_points_on_torus_d(
const Kernel &k,
int dim, std::size_t num_slices,
192 double radius_noise_percentage = 0.,
193 std::vector<typename Kernel::FT> current_point = std::vector<typename Kernel::FT>()) {
195 int point_size =
static_cast<int>(current_point.size());
196 if (point_size == 2 * dim) {
197 *out++ = k.construct_point_d_object()(point_size, current_point.begin(), current_point.end());
199 for (std::size_t slice_idx = 0; slice_idx < num_slices; ++slice_idx) {
200 double radius_noise_ratio = 1.;
201 if (radius_noise_percentage > 0.) {
202 radius_noise_ratio = rng.get_double(
203 (100. - radius_noise_percentage) / 100.,
204 (100. + radius_noise_percentage) / 100.);
206 std::vector<typename Kernel::FT> cp2 = current_point;
207 double alpha = 6.2832 * slice_idx / num_slices;
208 cp2.push_back(radius_noise_ratio * std::cos(alpha));
209 cp2.push_back(radius_noise_ratio * std::sin(alpha));
210 generate_uniform_points_on_torus_d(
211 k, dim, num_slices, out, radius_noise_percentage, cp2);
216 template <
typename Kernel>
217 std::vector<typename Kernel::Point_d> generate_points_on_torus_d(std::size_t num_points,
int dim,
bool uniform =
false,
218 double radius_noise_percentage = 0.) {
219 typedef typename Kernel::Point_d Point;
220 typedef typename Kernel::FT FT;
224 std::vector<Point> points;
225 points.reserve(num_points);
227 std::size_t num_slices = (std::size_t)std::pow(num_points, 1. / dim);
228 generate_uniform_points_on_torus_d(
229 k, dim, num_slices, std::back_inserter(points), radius_noise_percentage);
231 for (std::size_t i = 0; i < num_points;) {
232 double radius_noise_ratio = 1.;
233 if (radius_noise_percentage > 0.) {
234 radius_noise_ratio = rng.get_double(
235 (100. - radius_noise_percentage) / 100.,
236 (100. + radius_noise_percentage) / 100.);
238 std::vector<typename Kernel::FT> pt;
240 for (
int curdim = 0; curdim < dim; ++curdim) {
241 FT alpha = rng.get_double(0, 6.2832);
242 pt.push_back(radius_noise_ratio * std::cos(alpha));
243 pt.push_back(radius_noise_ratio * std::sin(alpha));
246 Point p = k.construct_point_d_object()(pt.begin(), pt.end());
254 template <
typename Kernel>
255 std::vector<typename Kernel::Point_d> generate_points_on_sphere_d(std::size_t num_points,
int dim,
double radius,
256 double radius_noise_percentage = 0.) {
257 typedef typename Kernel::Point_d Point;
260 CGAL::Random_points_on_sphere_d<Point> generator(dim, radius);
261 std::vector<Point> points;
262 points.reserve(num_points);
263 for (std::size_t i = 0; i < num_points;) {
264 Point p = *generator++;
265 if (radius_noise_percentage > 0.) {
266 double radius_noise_ratio = rng.get_double(
267 (100. - radius_noise_percentage) / 100.,
268 (100. + radius_noise_percentage) / 100.);
270 typename Kernel::Point_to_vector_d k_pt_to_vec =
271 k.point_to_vector_d_object();
272 typename Kernel::Vector_to_point_d k_vec_to_pt =
273 k.vector_to_point_d_object();
274 typename Kernel::Scaled_vector_d k_scaled_vec =
275 k.scaled_vector_d_object();
276 p = k_vec_to_pt(k_scaled_vec(k_pt_to_vec(p), radius_noise_ratio));
284 template <
typename Kernel>
285 std::vector<typename Kernel::Point_d> generate_points_in_ball_d(std::size_t num_points,
int dim,
double radius) {
286 typedef typename Kernel::Point_d Point;
289 CGAL::Random_points_in_ball_d<Point> generator(dim, radius);
290 std::vector<Point> points;
291 points.reserve(num_points);
292 for (std::size_t i = 0; i < num_points;) {
293 Point p = *generator++;
300 template <
typename Kernel>
301 std::vector<typename Kernel::Point_d> generate_points_in_cube_d(std::size_t num_points,
int dim,
double radius) {
302 typedef typename Kernel::Point_d Point;
305 CGAL::Random_points_in_cube_d<Point> generator(dim, radius);
306 std::vector<Point> points;
307 points.reserve(num_points);
308 for (std::size_t i = 0; i < num_points;) {
309 Point p = *generator++;
316 template <
typename Kernel>
317 std::vector<typename Kernel::Point_d> generate_points_on_two_spheres_d(std::size_t num_points,
int dim,
double radius,
318 double distance_between_centers,
319 double radius_noise_percentage = 0.) {
320 typedef typename Kernel::FT FT;
321 typedef typename Kernel::Point_d Point;
322 typedef typename Kernel::Vector_d Vector;
325 CGAL::Random_points_on_sphere_d<Point> generator(dim, radius);
326 std::vector<Point> points;
327 points.reserve(num_points);
329 std::vector<FT> t(dim, FT(0));
330 t[0] = distance_between_centers;
331 Vector c1_to_c2(t.begin(), t.end());
333 for (std::size_t i = 0; i < num_points;) {
334 Point p = *generator++;
335 if (radius_noise_percentage > 0.) {
336 double radius_noise_ratio = rng.get_double(
337 (100. - radius_noise_percentage) / 100.,
338 (100. + radius_noise_percentage) / 100.);
340 typename Kernel::Point_to_vector_d k_pt_to_vec =
341 k.point_to_vector_d_object();
342 typename Kernel::Vector_to_point_d k_vec_to_pt =
343 k.vector_to_point_d_object();
344 typename Kernel::Scaled_vector_d k_scaled_vec =
345 k.scaled_vector_d_object();
346 p = k_vec_to_pt(k_scaled_vec(k_pt_to_vec(p), radius_noise_ratio));
349 typename Kernel::Translated_point_d k_transl =
350 k.translated_point_d_object();
351 Point p2 = k_transl(p, c1_to_c2);
353 points.push_back(p2);
361 template <
typename Kernel>
362 std::vector<typename Kernel::Point_d> generate_points_on_3sphere_and_circle(std::size_t num_points,
363 double sphere_radius) {
364 typedef typename Kernel::FT FT;
365 typedef typename Kernel::Point_d Point;
368 CGAL::Random_points_on_sphere_d<Point> generator(3, sphere_radius);
369 std::vector<Point> points;
370 points.reserve(num_points);
372 typename Kernel::Compute_coordinate_d k_coord =
373 k.compute_coordinate_d_object();
374 for (std::size_t i = 0; i < num_points;) {
375 Point p_sphere = *generator++;
377 FT alpha = rng.get_double(0, 6.2832);
378 std::vector<FT> pt(5);
379 pt[0] = k_coord(p_sphere, 0);
380 pt[1] = k_coord(p_sphere, 1);
381 pt[2] = k_coord(p_sphere, 2);
382 pt[3] = std::cos(alpha);
383 pt[4] = std::sin(alpha);
384 Point p(pt.begin(), pt.end());
392 template <
typename Kernel>
393 std::vector<typename Kernel::Point_d> generate_points_on_klein_bottle_3D(std::size_t num_points,
double a,
double b,
394 bool uniform =
false) {
395 typedef typename Kernel::Point_d Point;
396 typedef typename Kernel::FT FT;
401 std::size_t num_lines = (std::size_t)sqrt(num_points);
403 std::vector<Point> points;
404 points.reserve(num_points);
405 for (std::size_t i = 0; i < num_points;) {
408 std::size_t k1 = i / num_lines;
409 std::size_t k2 = i % num_lines;
410 u = 6.2832 * k1 / num_lines;
411 v = 6.2832 * k2 / num_lines;
413 u = rng.get_double(0, 6.2832);
414 v = rng.get_double(0, 6.2832);
416 double tmp = cos(u / 2) * sin(v) - sin(u / 2) * sin(2. * v);
417 Point p = construct_point(k,
418 (a + b * tmp) * cos(u),
419 (a + b * tmp) * sin(u),
420 b * (sin(u / 2) * sin(v) + cos(u / 2) * sin(2. * v)));
428 template <
typename Kernel>
429 std::vector<typename Kernel::Point_d> generate_points_on_klein_bottle_4D(std::size_t num_points,
double a,
double b,
430 double noise = 0.,
bool uniform =
false) {
431 typedef typename Kernel::Point_d Point;
432 typedef typename Kernel::FT FT;
437 std::size_t num_lines = (std::size_t)sqrt(num_points);
439 std::vector<Point> points;
440 points.reserve(num_points);
441 for (std::size_t i = 0; i < num_points;) {
444 std::size_t k1 = i / num_lines;
445 std::size_t k2 = i % num_lines;
446 u = 6.2832 * k1 / num_lines;
447 v = 6.2832 * k2 / num_lines;
449 u = rng.get_double(0, 6.2832);
450 v = rng.get_double(0, 6.2832);
452 Point p = construct_point(k,
453 (a + b * cos(v)) * cos(u) + (noise == 0. ? 0. : rng.get_double(0, noise)),
454 (a + b * cos(v)) * sin(u) + (noise == 0. ? 0. : rng.get_double(0, noise)),
455 b * sin(v) * cos(u / 2) + (noise == 0. ? 0. : rng.get_double(0, noise)),
456 b * sin(v) * sin(u / 2) + (noise == 0. ? 0. : rng.get_double(0, noise)));
466 template <
typename Kernel>
467 std::vector<typename Kernel::Point_d>
468 generate_points_on_klein_bottle_variant_5D(
469 std::size_t num_points,
double a,
double b,
bool uniform =
false) {
470 typedef typename Kernel::Point_d Point;
471 typedef typename Kernel::FT FT;
476 std::size_t num_lines = (std::size_t)sqrt(num_points);
478 std::vector<Point> points;
479 points.reserve(num_points);
480 for (std::size_t i = 0; i < num_points;) {
483 std::size_t k1 = i / num_lines;
484 std::size_t k2 = i % num_lines;
485 u = 6.2832 * k1 / num_lines;
486 v = 6.2832 * k2 / num_lines;
488 u = rng.get_double(0, 6.2832);
489 v = rng.get_double(0, 6.2832);
491 FT x1 = (a + b * cos(v)) * cos(u);
492 FT x2 = (a + b * cos(v)) * sin(u);
493 FT x3 = b * sin(v) * cos(u / 2);
494 FT x4 = b * sin(v) * sin(u / 2);
495 FT x5 = x1 + x2 + x3 + x4;
497 Point p = construct_point(k, x1, x2, x3, x4, x5);
506 #endif // RANDOM_POINT_GENERATORS_H_ Definition: SimplicialComplexForAlpha.h:26