Implemented make_mesh_3 param :

initial_points_generator
This commit is contained in:
ange-clement 2023-10-13 14:51:43 +02:00
parent 8a4b492f1c
commit 3937c5df34
7 changed files with 462 additions and 230 deletions

View File

@ -167,6 +167,11 @@ if(TARGET CGAL::CGAL_ImageIO)
target_link_libraries(mesh_3D_gray_image_with_custom_initialization
PUBLIC CGAL::Eigen3_support)
create_single_source_cgal_program(
"mesh_3D_image_with_initial_points.cpp")
target_link_libraries(mesh_3D_image_with_initial_points
PUBLIC CGAL::Eigen3_support)
create_single_source_cgal_program("mesh_3D_image_variable_size.cpp")
target_link_libraries(mesh_3D_image_variable_size
PUBLIC CGAL::Eigen3_support)
@ -201,6 +206,7 @@ if(CGAL_ACTIVATE_CONCURRENT_MESH_3 AND TARGET CGAL::TBB_support)
mesh_3D_image_variable_size
mesh_3D_image_with_custom_initialization
mesh_3D_gray_image_with_custom_initialization
mesh_3D_image_with_initial_points
mesh_3D_image_with_features
mesh_3D_image_with_detection_of_features
mesh_3D_image_with_input_features

View File

@ -0,0 +1,58 @@
#include "random_labeled_image.h"
#include <CGAL/Exact_predicates_inexact_constructions_kernel.h>
#include <CGAL/Mesh_triangulation_3.h>
#include <CGAL/Mesh_complex_3_in_triangulation_3.h>
#include <CGAL/Mesh_criteria_3.h>
#include <CGAL/Mesh_3/Construct_initial_points_from_labeled_image.h>
#include <CGAL/Labeled_mesh_domain_3.h>
#include <CGAL/make_mesh_3.h>
#include <CGAL/Image_3.h>
#include <CGAL/SMDS_3/Dump_c3t3.h>
// Domain
typedef CGAL::Exact_predicates_inexact_constructions_kernel K;
typedef CGAL::Labeled_mesh_domain_3<K> Mesh_domain;
#ifdef CGAL_CONCURRENT_MESH_3
typedef CGAL::Parallel_tag Concurrency_tag;
#else
typedef CGAL::Sequential_tag Concurrency_tag;
#endif
// Triangulation
typedef CGAL::Mesh_triangulation_3<Mesh_domain,CGAL::Default,Concurrency_tag>::type Tr;
typedef CGAL::Mesh_complex_3_in_triangulation_3<Tr> C3t3;
// Criteria
typedef CGAL::Mesh_criteria_3<Tr> Mesh_criteria;
namespace params = CGAL::parameters;
int main()
{
/// [Create the image]
CGAL::Image_3 image = random_labeled_image();
/// [Create the image]
// Domain
Mesh_domain domain = Mesh_domain::create_labeled_image_mesh_domain(image);
// Mesh criteria
Mesh_criteria criteria(params::facet_angle(30).facet_size(3).facet_distance(1).
cell_radius_edge_ratio(3).cell_size(3));
C3t3 c3t3 = CGAL::make_mesh_3<C3t3>(domain, criteria,
params::initial_points_generator(Construct_initial_points_labeled_image(image))
);
/// [Meshing]
// Output
CGAL::dump_c3t3(c3t3, "out");
return 0;
}

View File

@ -0,0 +1,238 @@
// Copyright (c) 20XX,20XX GeometryFactory
// All rights reserved.
//
// This file is part of CGAL (www.cgal.org).
//
// $URL$
// $Id$
// SPDX-License-Identifier: GPL-3.0-or-later OR LicenseRef-Commercial
//
//
// Author(s) : Ange Clement
#ifndef CGAL_MESH_3_CONSTRUCT_INITIAL_POINTS_FROM_LABELED_IMAGE_H
#define CGAL_MESH_3_CONSTRUCT_INITIAL_POINTS_FROM_LABELED_IMAGE_H
#include <CGAL/license/Mesh_3.h>
#include <CGAL/Mesh_3/search_for_connected_components_in_labeled_image.h>
#include <CGAL/iterator.h>
#include <CGAL/point_generators_3.h>
#include <CGAL/Image_3.h>
template <typename Point>
struct Get_point
{
const double vx, vy, vz;
const double tx, ty, tz;
const std::size_t xdim, ydim, zdim;
Get_point(const CGAL::Image_3* image)
: vx(image->vx())
, vy(image->vy())
, vz(image->vz())
, tx(image->tx())
, ty(image->ty())
, tz(image->tz())
, xdim(image->xdim())
, ydim(image->ydim())
, zdim(image->zdim())
{}
Point operator()(const std::size_t i,
const std::size_t j,
const std::size_t k) const
{
double x = double(i) * vx + tx;
double y = double(j) * vy + ty;
double z = double(k) * vz + tz;
if (i == 0) x += 1. / 6. * vx;
else if (i == xdim - 1) x -= 1. / 6. * vx;
if (j == 0) y += 1. / 6. * vy;
else if (j == ydim - 1) y -= 1. / 6. * vy;
if (k == 0) z += 1. / 6. * vz;
else if (k == zdim - 1) z -= 1. / 6. * vz;
return Point(x, y, z);
}
};
struct Construct_initial_points_labeled_image
{
const CGAL::Image_3 & image;
Construct_initial_points_labeled_image(const CGAL::Image_3 & image_)
: image(image_)
{ }
template <typename OutputIterator, typename MeshDomain, typename C3t3>
OutputIterator operator()(OutputIterator pts, const MeshDomain& domain, const C3t3& c3t3, int n = 20) const
{
typedef typename MeshDomain::Subdomain Subdomain;
typedef typename MeshDomain::Point_3 Point_3;
typedef typename MeshDomain::Index Index;
typedef typename C3t3::Triangulation Tr;
typedef typename Tr::Geom_traits GT;
typedef typename GT::FT FT;
typedef typename Tr::Weighted_point Weighted_point;
typedef typename Tr::Segment Segment_3;
typedef typename Tr::Vertex_handle Vertex_handle;
typedef typename Tr::Cell_handle Cell_handle;
typedef typename GT::Vector_3 Vector_3;
const Tr& tr = c3t3.triangulation();
typename GT::Compare_weighted_squared_radius_3 cwsr =
tr.geom_traits().compare_weighted_squared_radius_3_object();
typename GT::Construct_point_3 cp =
tr.geom_traits().construct_point_3_object();
typename GT::Construct_weighted_point_3 cwp =
tr.geom_traits().construct_weighted_point_3_object();
const double max_v = (std::max)((std::max)(image.vx(),
image.vy()),
image.vz());
struct Seed {
std::size_t i, j, k;
std::size_t radius;
};
using Seeds = std::vector<Seed>;
Seeds seeds;
Get_point<Point_3> get_point(&image);
std::cout << "Searching for connected components..." << std::endl;
CGAL_IMAGE_IO_CASE(image.image(), search_for_connected_components_in_labeled_image(image,
std::back_inserter(seeds),
CGAL::Emptyset_iterator(),
CGAL::Identity<Word>(),
Word()));
std::cout << " " << seeds.size() << " components were found." << std::endl;
std::cout << "Construct initial points..." << std::endl;
for(const Seed seed : seeds)
{
const Point_3 seed_point = get_point(seed.i, seed.j, seed.k);
Cell_handle seed_cell = tr.locate(cwp(seed_point));
const Subdomain seed_label
= domain.is_in_domain_object()(seed_point);
const Subdomain seed_cell_label
= ( tr.dimension() < 3
|| seed_cell == Cell_handle()
|| tr.is_infinite(seed_cell))
? Subdomain() //seed_point is OUTSIDE_AFFINE_HULL
: domain.is_in_domain_object()(
seed_cell->weighted_circumcenter(tr.geom_traits()));
if ( seed_label != std::nullopt
&& seed_cell_label != std::nullopt
&& *seed_label == *seed_cell_label)
continue; //this means the connected component has already been initialized
const double radius = double(seed.radius + 1)* max_v;
CGAL::Random_points_on_sphere_3<Point_3> points_on_sphere_3(radius);
typename MeshDomain::Construct_intersection construct_intersection =
domain.construct_intersection_object();
std::vector<Vector_3> directions;
if(seed.radius < 2) {
// shoot in six directions
directions.push_back(Vector_3(-radius, 0, 0));
directions.push_back(Vector_3(+radius, 0, 0));
directions.push_back(Vector_3(0, -radius, 0));
directions.push_back(Vector_3(0, +radius, 0));
directions.push_back(Vector_3(0, 0, -radius));
directions.push_back(Vector_3(0, 0, +radius));
} else {
for(int i = 0; i < n; ++i)
{
// shoot n random directions
directions.push_back(*points_on_sphere_3++ - CGAL::ORIGIN);
}
}
for(const Vector_3& v : directions)
{
const Point_3 test = seed_point + v;
const Segment_3 test_segment = Segment_3(seed_point, test);
const typename MeshDomain::Intersection intersect =
construct_intersection(test_segment);
if (std::get<2>(intersect) != 0)
{
const Point_3& bpi = std::get<0>(intersect);
const Index index = std::get<1>(intersect);
Weighted_point pi = Weighted_point(bpi);
// This would cause trouble to optimizers
// check pi will not be hidden
typename Tr::Locate_type lt;
int li, lj;
Cell_handle pi_cell = tr.locate(pi, lt, li, lj);
if(lt != Tr::OUTSIDE_AFFINE_HULL) {
switch (tr.dimension())
{ //skip dimension 0
case 1:
if (tr.side_of_power_segment(pi_cell, pi, true) != CGAL::ON_BOUNDED_SIDE)
continue;
break;
case 2:
if (tr.side_of_power_circle(pi_cell, 3, pi, true) != CGAL::ON_BOUNDED_SIDE)
continue;
break;
case 3:
if (tr.side_of_power_sphere(pi_cell, pi, true) != CGAL::ON_BOUNDED_SIDE)
continue;
}
}
//check pi is not inside a protecting ball
std::vector<Vertex_handle> conflict_vertices;
if (tr.dimension() == 3)
{
tr.vertices_on_conflict_zone_boundary(pi, pi_cell
, std::back_inserter(conflict_vertices));
}
else
{
for (typename Tr::Finite_vertices_iterator vit = tr.finite_vertices_begin();
vit != tr.finite_vertices_end(); ++vit)
{
const Weighted_point& wp = tr.point(vit);
if (cwsr(wp, FT(0)) == CGAL::SMALLER) // 0 < wp's weight
conflict_vertices.push_back(vit);
}
}
bool pi_inside_protecting_sphere = false;
for(Vertex_handle cv : conflict_vertices)
{
if(tr.is_infinite(cv))
continue;
const Weighted_point& cv_wp = tr.point(cv);
if (cwsr(cv_wp, FT(0)) == CGAL::EQUAL) // 0 == wp's weight
continue;
// if the (squared) distance between bpi and cv is smaller or equal than cv's weight
if (cwsr(cv_wp, - tr.min_squared_distance(bpi, cp(cv_wp))) != CGAL::LARGER)
{
pi_inside_protecting_sphere = true;
break;
}
}
if (pi_inside_protecting_sphere)
continue;
*pts++ = std::make_pair(bpi, index);
}
}
}
return pts;
}
};
#endif // CGAL_MESH_3_CONSTRUCT_INITIAL_POINTS_FROM_LABELED_IMAGE_H

View File

@ -14,56 +14,10 @@
#define CGAL_MESH_3_INITIALIZE_TRIANGULATION_FROM_LABELED_IMAGE_H
#include <CGAL/license/Mesh_3.h>
#include <CGAL/Mesh_3/search_for_connected_components_in_labeled_image.h>
#include <CGAL/Distance_3/Point_3_Triangle_3.h>
#include <CGAL/Labeled_mesh_domain_3.h>
#include <CGAL/make_mesh_3.h>
#include <CGAL/enum.h>
#include <CGAL/iterator.h>
#include <CGAL/point_generators_3.h>
#include <CGAL/Image_3.h>
#include <CGAL/Mesh_3/Construct_initial_points_from_labeled_image.h>
#include <iostream>
#include <queue>
template <typename Point>
struct Get_point
{
const double vx, vy, vz;
const double tx, ty, tz;
const std::size_t xdim, ydim, zdim;
Get_point(const CGAL::Image_3* image)
: vx(image->vx())
, vy(image->vy())
, vz(image->vz())
, tx(image->tx())
, ty(image->ty())
, tz(image->tz())
, xdim(image->xdim())
, ydim(image->ydim())
, zdim(image->zdim())
{}
Point operator()(const std::size_t i,
const std::size_t j,
const std::size_t k) const
{
double x = double(i) * vx + tx;
double y = double(j) * vy + ty;
double z = double(k) * vz + tz;
if (i == 0) x += 1. / 6. * vx;
else if (i == xdim - 1) x -= 1. / 6. * vx;
if (j == 0) y += 1. / 6. * vy;
else if (j == ydim - 1) y -= 1. / 6. * vy;
if (k == 0) z += 1. / 6. * vz;
else if (k == zdim - 1) z -= 1. / 6. * vz;
return Point(x, y, z);
}
};
template<class C3T3, class MeshDomain, class MeshCriteria>
void init_tr_from_labeled_image_call_init_features(C3T3&,
const MeshDomain&,
@ -97,23 +51,17 @@ void initialize_triangulation_from_labeled_image(C3T3& c3t3,
{
typedef typename C3T3::Triangulation Tr;
typedef typename Tr::Geom_traits GT;
typedef typename GT::FT FT;
typedef typename Tr::Weighted_point Weighted_point;
typedef typename Tr::Bare_point Bare_point;
typedef typename Tr::Segment Segment_3;
typedef typename Tr::Vertex_handle Vertex_handle;
typedef typename Tr::Cell_handle Cell_handle;
typedef typename MeshDomain::Point_3 Point_3;
typedef typename MeshDomain::Index Index;
typedef typename GT::Vector_3 Vector_3;
typedef typename std::pair<Point_3, Index> ConstructedPoint;
typedef MeshDomain Mesh_domain;
Tr& tr = c3t3.triangulation();
typename GT::Compare_weighted_squared_radius_3 cwsr =
tr.geom_traits().compare_weighted_squared_radius_3_object();
typename GT::Construct_point_3 cp =
tr.geom_traits().construct_point_3_object();
typename GT::Construct_weighted_point_3 cwp =
tr.geom_traits().construct_weighted_point_3_object();
@ -123,159 +71,36 @@ void initialize_triangulation_from_labeled_image(C3T3& c3t3,
CGAL::internal::Has_features<Mesh_domain>());
}
const double max_v = (std::max)((std::max)(image.vx(),
image.vy()),
image.vz());
std::vector<ConstructedPoint> constructedPoints;
struct Seed {
std::size_t i, j, k;
std::size_t radius;
};
using Seeds = std::vector<Seed>;
using Subdomain = typename Mesh_domain::Subdomain;
Construct_initial_points_labeled_image construct(image);
construct(std::back_inserter(constructedPoints), domain, c3t3);
Seeds seeds;
Get_point<Bare_point> get_point(&image);
std::cout << "Searching for connected components..." << std::endl;
search_for_connected_components_in_labeled_image(image,
std::back_inserter(seeds),
CGAL::Emptyset_iterator(),
transform,
Image_word_type());
std::cout << " " << seeds.size() << " components were found." << std::endl;
std::cout << "Construct initial points..." << std::endl;
for(const Seed seed : seeds)
std::cout << " " << constructedPoints.size() << " constructed points" << std::endl;
for (const ConstructedPoint & constructedPoint : constructedPoints)
{
const Bare_point seed_point = get_point(seed.i, seed.j, seed.k);
Cell_handle seed_cell = tr.locate(cwp(seed_point));
const Point_3& point = constructedPoint.first;
const Index& index = constructedPoint.second;
const Subdomain seed_label
= domain.is_in_domain_object()(seed_point);
const Subdomain seed_cell_label
= ( tr.dimension() < 3
|| seed_cell == Cell_handle()
|| tr.is_infinite(seed_cell))
? Subdomain() //seed_point is OUTSIDE_AFFINE_HULL
: domain.is_in_domain_object()(
seed_cell->weighted_circumcenter(tr.geom_traits()));
Weighted_point pi = cwp(point);
if ( seed_label != std::nullopt
&& seed_cell_label != std::nullopt
&& *seed_label == *seed_cell_label)
continue; //this means the connected component has already been initialized
/// The following lines show how to insert initial points in the
/// `c3t3` object. [insert initial points]
Vertex_handle v = tr.insert(pi);
// `v` could be null if `pi` is hidden by other vertices of `tr`.
CGAL_assertion(v != Vertex_handle());
c3t3.set_dimension(v, 2); // by construction, points are on surface
c3t3.set_index(v, index);
/// [insert initial points]
}
const double radius = double(seed.radius + 1)* max_v;
CGAL::Random_points_on_sphere_3<Bare_point> points_on_sphere_3(radius);
typename Mesh_domain::Construct_intersection construct_intersection =
domain.construct_intersection_object();
std::vector<Vector_3> directions;
if(seed.radius < 2) {
// shoot in six directions
directions.push_back(Vector_3(-radius, 0, 0));
directions.push_back(Vector_3(+radius, 0, 0));
directions.push_back(Vector_3(0, -radius, 0));
directions.push_back(Vector_3(0, +radius, 0));
directions.push_back(Vector_3(0, 0, -radius));
directions.push_back(Vector_3(0, 0, +radius));
} else {
for(int i = 0; i < 20; ++i)
{
// shoot 20 random directions
directions.push_back(*points_on_sphere_3++ - CGAL::ORIGIN);
}
}
for(const Vector_3& v : directions)
{
const Bare_point test = seed_point + v;
const typename Mesh_domain::Intersection intersect =
construct_intersection(Segment_3(seed_point, test));
if (std::get<2>(intersect) != 0)
{
const Bare_point& bpi = std::get<0>(intersect);
Weighted_point pi = cwp(bpi);
// This would cause trouble to optimizers
// check pi will not be hidden
typename Tr::Locate_type lt;
int li, lj;
Cell_handle pi_cell = tr.locate(pi, lt, li, lj);
if(lt != Tr::OUTSIDE_AFFINE_HULL) {
switch (tr.dimension())
{ //skip dimension 0
case 1:
if (tr.side_of_power_segment(pi_cell, pi, true) != CGAL::ON_BOUNDED_SIDE)
continue;
break;
case 2:
if (tr.side_of_power_circle(pi_cell, 3, pi, true) != CGAL::ON_BOUNDED_SIDE)
continue;
break;
case 3:
if (tr.side_of_power_sphere(pi_cell, pi, true) != CGAL::ON_BOUNDED_SIDE)
continue;
}
}
//check pi is not inside a protecting ball
std::vector<Vertex_handle> conflict_vertices;
if (tr.dimension() == 3)
{
tr.vertices_on_conflict_zone_boundary(pi, pi_cell
, std::back_inserter(conflict_vertices));
}
else
{
for (typename Tr::Finite_vertices_iterator vit = tr.finite_vertices_begin();
vit != tr.finite_vertices_end(); ++vit)
{
const Weighted_point& wp = tr.point(vit);
if (cwsr(wp, FT(0)) == CGAL::SMALLER) // 0 < wp's weight
conflict_vertices.push_back(vit);
}
}
bool pi_inside_protecting_sphere = false;
for(Vertex_handle cv : conflict_vertices)
{
if(tr.is_infinite(cv))
continue;
const Weighted_point& cv_wp = tr.point(cv);
if (cwsr(cv_wp, FT(0)) == CGAL::EQUAL) // 0 == wp's weight
continue;
// if the (squared) distance between bpi and cv is smaller or equal than cv's weight
if (cwsr(cv_wp, - tr.min_squared_distance(bpi, cp(cv_wp))) != CGAL::LARGER)
{
pi_inside_protecting_sphere = true;
break;
}
}
if (pi_inside_protecting_sphere)
continue;
const typename Mesh_domain::Index index = std::get<1>(intersect);
/// The following lines show how to insert initial points in the
/// `c3t3` object. [insert initial points]
Vertex_handle v = tr.insert(pi);
// `v` could be null if `pi` is hidden by other vertices of `tr`.
CGAL_assertion(v != Vertex_handle());
c3t3.set_dimension(v, 2); // by construction, points are on surface
c3t3.set_index(v, index);
/// [insert initial points]
}
// else
// {
// std::cerr <<
// boost::format("Error. Segment (%1%, %2%) does not intersect the surface!\n")
// % it->first % test;
// }
}
if ( tr.dimension() != 3 )
{
std::cout << " not enough points: triangulation.dimension() == "
<< tr.dimension() << std::endl;
CGAL::Mesh_3::internal::init_c3t3(c3t3, domain, criteria, 20);
std::cout << " -> " << tr.number_of_vertices() << " initial points." << std::endl;
}
std::cout << " " << tr.number_of_vertices() << " initial points." << std::endl;
if ( c3t3.triangulation().dimension() != 3 )

View File

@ -38,10 +38,10 @@ namespace CGAL {
namespace Mesh_3 {
namespace internal {
template < typename C3T3, typename MeshDomain, typename MeshCriteria >
template < typename C3T3, typename MeshDomain, typename MeshCriteria, typename InitialPointsGenerator = Null_functor >
void
init_c3t3(C3T3& c3t3, const MeshDomain& domain, const MeshCriteria&,
const int nb_initial_points)
const int nb_initial_points, InitialPointsGenerator& generator = Null_functor())
{
typedef typename MeshDomain::Point_3 Point_3;
typedef typename MeshDomain::Index Index;
@ -49,13 +49,23 @@ init_c3t3(C3T3& c3t3, const MeshDomain& domain, const MeshCriteria&,
typedef typename Initial_points_vector::iterator Ipv_iterator;
typedef typename C3T3::Vertex_handle Vertex_handle;
typedef typename
CGAL::parameters::internal::Initial_points_generator_generator<
std::back_insert_iterator<Initial_points_vector>,
MeshDomain,
C3T3,
InitialPointsGenerator>
Initial_points_generator_generator;
// Mesh initialization : get some points and add them to the mesh
Initial_points_vector initial_points;
auto& generator_wrapped =
Initial_points_generator_generator()(generator);
if (nb_initial_points > -1)
domain.construct_initial_points_object()(std::back_inserter(initial_points),
generator_wrapped(std::back_inserter(initial_points), domain, c3t3,
nb_initial_points);
else //use default number of points
domain.construct_initial_points_object()(std::back_inserter(initial_points));
generator_wrapped(std::back_inserter(initial_points), domain, c3t3);
typename C3T3::Triangulation::Geom_traits::Construct_weighted_point_3 cwp =
c3t3.triangulation().geom_traits().construct_weighted_point_3_object();
@ -159,20 +169,22 @@ template < typename C3T3,
typename MeshDomain,
typename MeshCriteria,
bool MeshDomainHasHasFeatures,
typename HasFeatures = int>
typename HasFeatures = int,
typename InitialPointsGenerator = Null_functor>
struct C3t3_initializer { };
// Partial specialization of C3t3_initializer
// Handle cases where MeshDomain::Has_features is not a valid type
template < typename C3T3, typename MD, typename MC, typename HasFeatures >
struct C3t3_initializer < C3T3, MD, MC, false, HasFeatures >
template < typename C3T3, typename MD, typename MC, typename HasFeatures, typename InitialPointsGenerator >
struct C3t3_initializer < C3T3, MD, MC, false, HasFeatures, InitialPointsGenerator >
{
typedef parameters::internal::Mesh_3_options Mesh_3_options;
void operator()(C3T3& c3t3,
const MD& domain,
const MC& criteria,
bool with_features,
Mesh_3_options mesh_options = Mesh_3_options())
Mesh_3_options mesh_options = Mesh_3_options(),
InitialPointsGenerator& generator = Null_functor())
{
if ( with_features )
{
@ -181,23 +193,24 @@ struct C3t3_initializer < C3T3, MD, MC, false, HasFeatures >
}
init_c3t3(c3t3,domain,criteria,
mesh_options.number_of_initial_points);
mesh_options.number_of_initial_points, generator);
}
};
// Partial specialization of C3t3_initializer
// Handles cases where MeshDomain::Has_features is a valid type
template < typename C3T3, typename MD, typename MC, typename HasFeatures >
struct C3t3_initializer < C3T3, MD, MC, true, HasFeatures >
template < typename C3T3, typename MD, typename MC, typename HasFeatures, typename InitialPointsGenerator >
struct C3t3_initializer < C3T3, MD, MC, true, HasFeatures, InitialPointsGenerator >
{
typedef parameters::internal::Mesh_3_options Mesh_3_options;
void operator()(C3T3& c3t3,
const MD& domain,
const MC& criteria,
bool with_features,
Mesh_3_options mesh_options = Mesh_3_options())
Mesh_3_options mesh_options = Mesh_3_options(),
InitialPointsGenerator& generator = Null_functor())
{
C3t3_initializer < C3T3, MD, MC, true, typename MD::Has_features >()
C3t3_initializer < C3T3, MD, MC, true, typename MD::Has_features, InitialPointsGenerator >()
(c3t3,domain,criteria,with_features,mesh_options);
}
};
@ -205,8 +218,8 @@ struct C3t3_initializer < C3T3, MD, MC, true, HasFeatures >
// Partial specialization of C3t3_initializer
// Handles cases where MeshDomain::Has_features is a valid type and is defined
// to CGAL::Tag_true
template < typename C3T3, typename MD, typename MC >
struct C3t3_initializer < C3T3, MD, MC, true, CGAL::Tag_true >
template < typename C3T3, typename MD, typename MC, typename InitialPointsGenerator >
struct C3t3_initializer < C3T3, MD, MC, true, CGAL::Tag_true, InitialPointsGenerator >
: public C3t3_initializer_base < C3T3, MD, MC >
{
virtual ~C3t3_initializer() { }
@ -216,7 +229,8 @@ struct C3t3_initializer < C3T3, MD, MC, true, CGAL::Tag_true >
const MD& domain,
const MC& criteria,
bool with_features,
Mesh_3_options mesh_options = Mesh_3_options())
Mesh_3_options mesh_options = Mesh_3_options(),
InitialPointsGenerator& generator = Null_functor())
{
if ( with_features ) {
this->initialize_features(c3t3, domain, criteria,mesh_options);
@ -242,26 +256,27 @@ struct C3t3_initializer < C3T3, MD, MC, true, CGAL::Tag_true >
}
if(need_more_init) {
init_c3t3(c3t3, domain, criteria,
mesh_options.number_of_initial_points);
mesh_options.number_of_initial_points, generator);
}
}
else { init_c3t3(c3t3,domain,criteria,
mesh_options.number_of_initial_points); }
mesh_options.number_of_initial_points, generator); }
}
};
// Partial specialization of C3t3_initializer
// Handles cases where MeshDomain::Has_features is a valid type and is defined
// to CGAL::Tag_false
template < typename C3T3, typename MD, typename MC >
struct C3t3_initializer < C3T3, MD, MC, true, CGAL::Tag_false >
template < typename C3T3, typename MD, typename MC, typename InitialPointsGenerator >
struct C3t3_initializer < C3T3, MD, MC, true, CGAL::Tag_false, InitialPointsGenerator >
{
typedef parameters::internal::Mesh_3_options Mesh_3_options;
void operator()(C3T3& c3t3,
const MD& domain,
const MC& criteria,
bool with_features,
Mesh_3_options mesh_options = Mesh_3_options())
Mesh_3_options mesh_options = Mesh_3_options(),
InitialPointsGenerator& generator = Null_functor())
{
if ( with_features )
{
@ -270,7 +285,7 @@ struct C3t3_initializer < C3T3, MD, MC, true, CGAL::Tag_false >
}
init_c3t3(c3t3,domain,criteria,
mesh_options.number_of_initial_points);
mesh_options.number_of_initial_points, generator);
}
};
@ -401,6 +416,31 @@ struct C3t3_initializer < C3T3, MD, MC, true, CGAL::Tag_false >
* </UL>}
* \cgalParamDefault{`parameters::exude()`}
* \cgalParamSectionEnd
* \cgalParamSectionBegin{Mesh initialisation}
* \cgalParamDescription{an `InitialPointsGenerator` can optionally be supplied before the meshing process.
* It must folow the `InitialPointsGenerator` concept and is specified with the param:
* <UL>
* <LI> `parameters::initial_points_generator(generator)`
* </UL>
*
* The `InitialPointsGenerator` concept is a function object to construct
* a set of initial points on the surface of the domain. Provides the
* following operators:
*
* `template <typename OutputIterator, typename MeshDomain, typename C3t3>`
* <br>
* `OutputIterator operator()(OutputIterator pts, const MeshDomain& domain, const C3t3& c3t3)`
*
* `template <typename OutputIterator, typename MeshDomain, typename C3t3>`
* <br>
* `OutputIterator operator()(OutputIterator pts, const MeshDomain& domain, const C3t3& c3t3, int n)`
*
* Those two operators output a set of (`n`) surface points to the
* output iterator `pts`, as objects of type `std::pair<Point_3,
* %Index>`. If `n` is not given, the functor must provide enough
* points to initialize the mesh generation process.}
* \cgalParamDefault{`parameters::initial_points_generator(generator)`}
* \cgalParamSectionEnd
* \cgalNamedParamsEnd
*
* Note that regardless of which optimization processes are activated,
@ -436,11 +476,13 @@ C3T3 make_mesh_3(const MeshDomain& domain, const MeshCriteria& criteria, const C
parameters::internal::Features_options features_param = choose_parameter(get_parameter(np, internal_np::features_options_param), parameters::features(domain).v);
parameters::internal::Mesh_3_options mesh_options_param = choose_parameter(get_parameter(np, internal_np::mesh_param), parameters::internal::Mesh_3_options());
parameters::internal::Manifold_options manifold_options_param = choose_parameter(get_parameter(np, internal_np::manifold_param), parameters::internal::Manifold_options());
auto initial_points_generator_param = choose_parameter(get_parameter(np, internal_np::initial_points_generator_param), Null_functor());
make_mesh_3_impl(c3t3, domain, criteria,
exude_param, perturb_param, odt_param, lloyd_param,
features_param.features(), mesh_options_param,
manifold_options_param);
manifold_options_param,
initial_points_generator_param);
return c3t3;
}
@ -469,7 +511,7 @@ C3T3 make_mesh_3(const MeshDomain& domain, const MeshCriteria& criteria,
*
* @return The mesh as a C3T3 object
*/
template<class C3T3, class MeshDomain, class MeshCriteria>
template<class C3T3, class MeshDomain, class MeshCriteria, class InitialPointsGenerator = Null_functor>
void make_mesh_3_impl(C3T3& c3t3,
const MeshDomain& domain,
const MeshCriteria& criteria,
@ -481,7 +523,8 @@ void make_mesh_3_impl(C3T3& c3t3,
const parameters::internal::Mesh_3_options&
mesh_options = parameters::internal::Mesh_3_options(),
const parameters::internal::Manifold_options&
manifold_options = parameters::internal::Manifold_options())
manifold_options = parameters::internal::Manifold_options(),
InitialPointsGenerator& generator = Null_functor())
{
#ifdef CGAL_MESH_3_INITIAL_POINTS_NO_RANDOM_SHOOTING
CGAL::get_default_random() = CGAL::Random(0);
@ -492,11 +535,14 @@ void make_mesh_3_impl(C3T3& c3t3,
C3T3,
MeshDomain,
MeshCriteria,
::CGAL::internal::has_Has_features<MeshDomain>::value > () (c3t3,
::CGAL::internal::has_Has_features<MeshDomain>::value,
int, /*Type of MeshDomain::Has_features not determined*/
InitialPointsGenerator >() (c3t3,
domain,
criteria,
with_features,
mesh_options);
mesh_options,
generator);
CGAL_assertion( c3t3.triangulation().dimension() >= 2 );

View File

@ -11,6 +11,8 @@
#ifndef CGAL_MESH_OPTION_CLASSES_H
#define CGAL_MESH_OPTION_CLASSES_H
#include <functional>
#include <CGAL/STL_Extension/internal/Has_features.h>
namespace CGAL {
@ -165,6 +167,43 @@ private:
bool b_;
};
// Initial points generator
template <typename OutputIterator, typename MeshDomain, typename C3t3>
struct Initial_points_generator_wrapper
{
template <typename Initial_points_generator>
Initial_points_generator_wrapper(Initial_points_generator& generator)
: initial_points_generator_default_(generator)
, initial_points_generator_(generator)
{ }
OutputIterator operator()(OutputIterator pts, const MeshDomain& domain, const C3t3& c3t3)
{
return initial_points_generator_default_(pts, domain, c3t3);
}
OutputIterator operator()(OutputIterator pts, const MeshDomain& domain, const C3t3& c3t3, int n)
{
return initial_points_generator_(pts, domain, c3t3, n);
}
private:
std::function<OutputIterator(OutputIterator,MeshDomain,C3t3)> initial_points_generator_default_;
std::function<OutputIterator(OutputIterator,MeshDomain,C3t3,int)> initial_points_generator_;
};
template <typename OutputIterator, typename MeshDomain, typename C3t3>
struct Initial_points_generator_default
{
OutputIterator operator()(OutputIterator pts, const MeshDomain& domain, const C3t3& c3t3)
{
return domain.construct_initial_points_object()(pts);
}
OutputIterator operator()(OutputIterator pts, const MeshDomain& domain, const C3t3& c3t3, int n)
{
return domain.construct_initial_points_object()(pts, n);
}
};
// -----------------------------------
// Features generator
// -----------------------------------
@ -207,6 +246,25 @@ struct Domain_features_generator< MeshDomain, true >
}
};
// struct Initial_points_generator_generator
template <typename OutputIterator, typename MeshDomain, typename C3t3, typename Initial_points_generator>
struct Initial_points_generator_generator
{
auto operator()(Initial_points_generator& generator)
{
return Initial_points_generator_wrapper<OutputIterator, MeshDomain, C3t3>(generator);
}
};
template <typename OutputIterator, typename MeshDomain, typename C3t3>
struct Initial_points_generator_generator<OutputIterator, MeshDomain, C3t3, Null_functor>
{
auto operator()(Null_functor&)
{
return Initial_points_generator_default<OutputIterator, MeshDomain, C3t3>();
}
};
} // end namespace internal

View File

@ -323,6 +323,7 @@ CGAL_add_named_parameter_with_compatibility(do_reset_c3t3_t, do_reset_c3t3, do_r
CGAL_add_named_parameter_with_compatibility(mesh_param_t, mesh_param, mesh_options)
CGAL_add_named_parameter_with_compatibility(manifold_param_t, manifold_param, manifold_option)
CGAL_add_named_parameter_with_compatibility(features_option_param_t,features_options_param,features_options)
CGAL_add_named_parameter_with_compatibility(initial_points_generator_param_t,initial_points_generator_param,initial_points_generator)
CGAL_add_named_parameter_with_compatibility_cref_only(image_3_param_t, image_3_param, image)
CGAL_add_named_parameter_with_compatibility(iso_value_param_t, iso_value_param, iso_value)