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 =
194 std::vector<typename Kernel::FT>()) {
196 int point_size =
static_cast<int>(current_point.size());
197 if (point_size == 2 * dim) {
198 *out++ = k.construct_point_d_object()(point_size, current_point.begin(), current_point.end());
200 for (std::size_t slice_idx = 0; slice_idx < num_slices; ++slice_idx) {
201 double radius_noise_ratio = 1.;
202 if (radius_noise_percentage > 0.) {
203 radius_noise_ratio = rng.get_double(
204 (100. - radius_noise_percentage) / 100.,
205 (100. + radius_noise_percentage) / 100.);
207 std::vector<typename Kernel::FT> cp2 = current_point;
208 double alpha = 6.2832 * slice_idx / num_slices;
209 cp2.push_back(radius_noise_ratio * std::cos(alpha));
210 cp2.push_back(radius_noise_ratio * std::sin(alpha));
211 generate_uniform_points_on_torus_d(
212 k, dim, num_slices, out, radius_noise_percentage, cp2);
217 template <
typename Kernel>
218 std::vector<typename Kernel::Point_d> generate_points_on_torus_d(std::size_t num_points,
int dim,
bool uniform =
false,
219 double radius_noise_percentage = 0.) {
220 typedef typename Kernel::Point_d Point;
221 typedef typename Kernel::FT FT;
225 std::vector<Point> points;
226 points.reserve(num_points);
228 std::size_t num_slices = (std::size_t)std::pow(num_points, 1. / dim);
229 generate_uniform_points_on_torus_d(
230 k, dim, num_slices, std::back_inserter(points), radius_noise_percentage);
232 for (std::size_t i = 0; i < num_points;) {
233 double radius_noise_ratio = 1.;
234 if (radius_noise_percentage > 0.) {
235 radius_noise_ratio = rng.get_double(
236 (100. - radius_noise_percentage) / 100.,
237 (100. + radius_noise_percentage) / 100.);
239 std::vector<typename Kernel::FT> pt;
241 for (
int curdim = 0; curdim < dim; ++curdim) {
242 FT alpha = rng.get_double(0, 6.2832);
243 pt.push_back(radius_noise_ratio * std::cos(alpha));
244 pt.push_back(radius_noise_ratio * std::sin(alpha));
247 Point p = k.construct_point_d_object()(pt.begin(), pt.end());
255 template <
typename Kernel>
256 std::vector<typename Kernel::Point_d> generate_points_on_sphere_d(std::size_t num_points,
int dim,
double radius,
257 double radius_noise_percentage = 0.) {
258 typedef typename Kernel::Point_d Point;
261 CGAL::Random_points_on_sphere_d<Point> generator(dim, radius);
262 std::vector<Point> points;
263 points.reserve(num_points);
264 for (std::size_t i = 0; i < num_points;) {
265 Point p = *generator++;
266 if (radius_noise_percentage > 0.) {
267 double radius_noise_ratio = rng.get_double(
268 (100. - radius_noise_percentage) / 100.,
269 (100. + radius_noise_percentage) / 100.);
271 typename Kernel::Point_to_vector_d k_pt_to_vec =
272 k.point_to_vector_d_object();
273 typename Kernel::Vector_to_point_d k_vec_to_pt =
274 k.vector_to_point_d_object();
275 typename Kernel::Scaled_vector_d k_scaled_vec =
276 k.scaled_vector_d_object();
277 p = k_vec_to_pt(k_scaled_vec(k_pt_to_vec(p), radius_noise_ratio));
285 template <
typename Kernel>
286 std::vector<typename Kernel::Point_d> generate_points_in_ball_d(std::size_t num_points,
int dim,
double radius) {
287 typedef typename Kernel::Point_d Point;
290 CGAL::Random_points_in_ball_d<Point> generator(dim, radius);
291 std::vector<Point> points;
292 points.reserve(num_points);
293 for (std::size_t i = 0; i < num_points;) {
294 Point p = *generator++;
301 template <
typename Kernel>
302 std::vector<typename Kernel::Point_d> generate_points_in_cube_d(std::size_t num_points,
int dim,
double radius) {
303 typedef typename Kernel::Point_d Point;
306 CGAL::Random_points_in_cube_d<Point> generator(dim, radius);
307 std::vector<Point> points;
308 points.reserve(num_points);
309 for (std::size_t i = 0; i < num_points;) {
310 Point p = *generator++;
317 template <
typename Kernel>
318 std::vector<typename Kernel::Point_d> generate_points_on_two_spheres_d(std::size_t num_points,
int dim,
double radius,
319 double distance_between_centers,
320 double radius_noise_percentage = 0.) {
321 typedef typename Kernel::FT FT;
322 typedef typename Kernel::Point_d Point;
323 typedef typename Kernel::Vector_d Vector;
326 CGAL::Random_points_on_sphere_d<Point> generator(dim, radius);
327 std::vector<Point> points;
328 points.reserve(num_points);
330 std::vector<FT> t(dim, FT(0));
331 t[0] = distance_between_centers;
332 Vector c1_to_c2(t.begin(), t.end());
334 for (std::size_t i = 0; i < num_points;) {
335 Point p = *generator++;
336 if (radius_noise_percentage > 0.) {
337 double radius_noise_ratio = rng.get_double(
338 (100. - radius_noise_percentage) / 100.,
339 (100. + radius_noise_percentage) / 100.);
341 typename Kernel::Point_to_vector_d k_pt_to_vec =
342 k.point_to_vector_d_object();
343 typename Kernel::Vector_to_point_d k_vec_to_pt =
344 k.vector_to_point_d_object();
345 typename Kernel::Scaled_vector_d k_scaled_vec =
346 k.scaled_vector_d_object();
347 p = k_vec_to_pt(k_scaled_vec(k_pt_to_vec(p), radius_noise_ratio));
350 typename Kernel::Translated_point_d k_transl =
351 k.translated_point_d_object();
352 Point p2 = k_transl(p, c1_to_c2);
354 points.push_back(p2);
362 template <
typename Kernel>
363 std::vector<typename Kernel::Point_d> generate_points_on_3sphere_and_circle(std::size_t num_points,
364 double sphere_radius) {
365 typedef typename Kernel::FT FT;
366 typedef typename Kernel::Point_d Point;
369 CGAL::Random_points_on_sphere_d<Point> generator(3, sphere_radius);
370 std::vector<Point> points;
371 points.reserve(num_points);
373 typename Kernel::Compute_coordinate_d k_coord =
374 k.compute_coordinate_d_object();
375 for (std::size_t i = 0; i < num_points;) {
376 Point p_sphere = *generator++;
378 FT alpha = rng.get_double(0, 6.2832);
379 std::vector<FT> pt(5);
380 pt[0] = k_coord(p_sphere, 0);
381 pt[1] = k_coord(p_sphere, 1);
382 pt[2] = k_coord(p_sphere, 2);
383 pt[3] = std::cos(alpha);
384 pt[4] = std::sin(alpha);
385 Point p(pt.begin(), pt.end());
393 template <
typename Kernel>
394 std::vector<typename Kernel::Point_d> generate_points_on_klein_bottle_3D(std::size_t num_points,
double a,
double b,
395 bool uniform =
false) {
396 typedef typename Kernel::Point_d Point;
397 typedef typename Kernel::FT FT;
402 std::size_t num_lines = (std::size_t)sqrt(num_points);
404 std::vector<Point> points;
405 points.reserve(num_points);
406 for (std::size_t i = 0; i < num_points;) {
409 std::size_t k1 = i / num_lines;
410 std::size_t k2 = i % num_lines;
411 u = 6.2832 * k1 / num_lines;
412 v = 6.2832 * k2 / num_lines;
414 u = rng.get_double(0, 6.2832);
415 v = rng.get_double(0, 6.2832);
417 double tmp = cos(u / 2) * sin(v) - sin(u / 2) * sin(2. * v);
418 Point p = construct_point(k,
419 (a + b * tmp) * cos(u),
420 (a + b * tmp) * sin(u),
421 b * (sin(u / 2) * sin(v) + cos(u / 2) * sin(2. * v)));
429 template <
typename Kernel>
430 std::vector<typename Kernel::Point_d> generate_points_on_klein_bottle_4D(std::size_t num_points,
double a,
double b,
431 double noise = 0.,
bool uniform =
false) {
432 typedef typename Kernel::Point_d Point;
433 typedef typename Kernel::FT FT;
438 std::size_t num_lines = (std::size_t)sqrt(num_points);
440 std::vector<Point> points;
441 points.reserve(num_points);
442 for (std::size_t i = 0; i < num_points;) {
445 std::size_t k1 = i / num_lines;
446 std::size_t k2 = i % num_lines;
447 u = 6.2832 * k1 / num_lines;
448 v = 6.2832 * k2 / num_lines;
450 u = rng.get_double(0, 6.2832);
451 v = rng.get_double(0, 6.2832);
453 Point p = construct_point(k,
454 (a + b * cos(v)) * cos(u) + (noise == 0. ? 0. : rng.get_double(0, noise)),
455 (a + b * cos(v)) * sin(u) + (noise == 0. ? 0. : rng.get_double(0, noise)),
456 b * sin(v) * cos(u / 2) + (noise == 0. ? 0. : rng.get_double(0, noise)),
457 b * sin(v) * sin(u / 2) + (noise == 0. ? 0. : rng.get_double(0, noise)));
467 template <
typename Kernel>
468 std::vector<typename Kernel::Point_d>
469 generate_points_on_klein_bottle_variant_5D(
470 std::size_t num_points,
double a,
double b,
bool uniform =
false) {
471 typedef typename Kernel::Point_d Point;
472 typedef typename Kernel::FT FT;
477 std::size_t num_lines = (std::size_t)sqrt(num_points);
479 std::vector<Point> points;
480 points.reserve(num_points);
481 for (std::size_t i = 0; i < num_points;) {
484 std::size_t k1 = i / num_lines;
485 std::size_t k2 = i % num_lines;
486 u = 6.2832 * k1 / num_lines;
487 v = 6.2832 * k2 / num_lines;
489 u = rng.get_double(0, 6.2832);
490 v = rng.get_double(0, 6.2832);
492 FT x1 = (a + b * cos(v)) * cos(u);
493 FT x2 = (a + b * cos(v)) * sin(u);
494 FT x3 = b * sin(v) * cos(u / 2);
495 FT x4 = b * sin(v) * sin(u / 2);
496 FT x5 = x1 + x2 + x3 + x4;
498 Point p = construct_point(k, x1, x2, x3, x4, x5);
507 #endif // RANDOM_POINT_GENERATORS_H_ Definition: SimplicialComplexForAlpha.h:26