CentroidalVoronoiTesselation is the main component of the remeshing algorithm.
More...
#include <geogram/voronoi/CVT.h>
|
| CentroidalVoronoiTesselation (Mesh *mesh, coord_index_t dimension=0, const std::string &delaunay="default") |
| Constructs a new CentroidalVoronoiTesselation. More...
|
|
| CentroidalVoronoiTesselation (Mesh *mesh, const vector< vec3 > &R3_embedding, coord_index_t dimension=0, const std::string &delaunay="default") |
| Constructs a new CentroidalVoronoiTesselation. More...
|
|
virtual | ~CentroidalVoronoiTesselation () |
| Destructor.
|
|
bool | compute_initial_sampling (index_t nb_samples, bool verbose=false) |
| Computes a random initial sampling of the surface in nD. More...
|
|
void | set_points (index_t nb_points, const double *points) |
| Initializes the points with a user-specified vector. More...
|
|
void | resize_points (index_t nb_points) |
| Changes the number of points. More...
|
|
virtual void | Lloyd_iterations (index_t nb_iter) |
| Relaxes the points with Lloyd's algorithm. More...
|
|
virtual void | Newton_iterations (index_t nb_iter, index_t m=7) |
| Relaxes the points with Newton-Lloyd's algorithm. More...
|
|
void | compute_surface (Mesh *mesh, bool multinerve=true) |
| Computes the surfacic mesh (using the current points). More...
|
|
void | compute_volume (Mesh *mesh) |
| Computes the volumetric mesh (using the current points). More...
|
|
void | set_show_iterations (bool x) |
| Specifies whether a progress bar should be used. More...
|
|
void | set_use_RVC_centroids (bool x) |
| Specifies whether centroids of Voronoi cells should be used. More...
|
|
void | set_constrained_cvt (bool x) |
| Specifies whether constrained mode should be used. More...
|
|
Mesh * | mesh () |
|
Delaunay * | delaunay () |
|
RestrictedVoronoiDiagram * | RVD () |
|
void | set_facets_range (index_t facets_begin, index_t facets_end) |
| Restricts computation to a part of the input mesh. More...
|
|
void | make_current () |
| Makes this CentroidalVoronoiTesselation the current one. More...
|
|
void | done_current () |
| Resets the current CentroidalVoronoiTesselation to nullptr. More...
|
|
void | set_progress_logger (ProgressTask *progress) |
| Sets a client for the progress bars. More...
|
|
coord_index_t | dimension () const |
| Gets the dimension of the points. More...
|
|
index_t | nb_points () const |
| Gets the number of points to be optimized.
|
|
const vec3 & | R3_embedding (index_t p) const |
| Gets the representation of a point in R3. More...
|
|
double * | embedding (index_t p) |
| Returns the representation of a point in embedding space. More...
|
|
bool | volumetric () const |
| Tests whether volumetric mode is used.
|
|
void | set_volumetric (bool x) |
| Sets volumetric mode. More...
|
|
bool | point_is_locked (index_t i) const |
| Tests whether a point is locked. More...
|
|
void | lock_point (index_t i) |
| Locks a point. More...
|
|
void | unlock_point (index_t i) |
| Unlocks a point. More...
|
|
void | unlock_all_points () |
| Unlocks all the points. More...
|
|
CentroidalVoronoiTesselation is the main component of the remeshing algorithm.
It evenly distributes points over a surface embedded in Rn, where n can be of arbitrary dimension. The geometrical computations are done by RestrictedVoronoiDiagram, and the numerical optimization by Optimizer.
Definition at line 69 of file CVT.h.
◆ CentroidalVoronoiTesselation() [1/2]
GEO::CentroidalVoronoiTesselation::CentroidalVoronoiTesselation |
( |
Mesh * |
mesh, |
|
|
coord_index_t |
dimension = 0 , |
|
|
const std::string & |
delaunay = "default" |
|
) |
| |
Constructs a new CentroidalVoronoiTesselation.
This constructor should be used when the first three coordinates of the mesh are x,y,z.
- Parameters
-
[in] | mesh | a pointer to the input mesh |
[in] | dimension | If set, uses only the dimension first coordinates in mesh, else dimension is determined by mesh->dimension(). |
[in] | delaunay | factory name of the implementation of Delaunay triangulation. Default uses ANN and radius of security. |
◆ CentroidalVoronoiTesselation() [2/2]
GEO::CentroidalVoronoiTesselation::CentroidalVoronoiTesselation |
( |
Mesh * |
mesh, |
|
|
const vector< vec3 > & |
R3_embedding, |
|
|
coord_index_t |
dimension = 0 , |
|
|
const std::string & |
delaunay = "default" |
|
) |
| |
Constructs a new CentroidalVoronoiTesselation.
This constructor should be used when the coordinates of the mesh are not related with R3.
- Parameters
-
[in] | mesh | a pointer to the input mesh |
[in] | R3_embedding | (dimension = mesh->nb_vertices()): coordinates of the mesh vertices in R3. Ignored if size is zero. |
[in] | dimension | If set, uses only the dimension first coordinates in mesh, else dimension is determined by mesh->dimension(). |
[in] | delaunay | factory name of the implementation of Delaunay triangulation. delaunay="default" uses ANN and radius of security. |
◆ compute_initial_sampling()
bool GEO::CentroidalVoronoiTesselation::compute_initial_sampling |
( |
index_t |
nb_samples, |
|
|
bool |
verbose = false |
|
) |
| |
Computes a random initial sampling of the surface in nD.
This initial sampling (of low quality/regularity) needs to be further optimized (using Lloyd_iterations() and Newton_iterations()).
- Parameters
-
[in] | nb_samples | number of points to generate in the sampling |
[in] | verbose | if set, display message |
◆ compute_R3_embedding()
void GEO::CentroidalVoronoiTesselation::compute_R3_embedding |
( |
| ) |
|
|
protected |
Computes the 3d representation of the Nd points.
It projects the points onto the Nd surface, then recovers the 3d coordinates by barycentric interpolation.
◆ compute_surface()
void GEO::CentroidalVoronoiTesselation::compute_surface |
( |
Mesh * |
mesh, |
|
|
bool |
multinerve = true |
|
) |
| |
Computes the surfacic mesh (using the current points).
- Parameters
-
[out] | mesh | the computed surface |
[in] | multinerve | If set, does topology control (uses the dual of the connected components of the RVD). |
◆ compute_volume()
void GEO::CentroidalVoronoiTesselation::compute_volume |
( |
Mesh * |
mesh | ) |
|
Computes the volumetric mesh (using the current points).
- Parameters
-
[out] | mesh | the computed volumetric mesh |
- Precondition
- volumetric()
◆ constrain_points()
void GEO::CentroidalVoronoiTesselation::constrain_points |
( |
double * |
g | ) |
const |
|
protected |
Constrains the locked points.
Zeroes the gradient relative to the components of locked points.
- Parameters
-
[in,out] | g | gradient of the objective function |
◆ delaunay()
Delaunay* GEO::CentroidalVoronoiTesselation::delaunay |
( |
| ) |
|
|
inline |
◆ dimension()
coord_index_t GEO::CentroidalVoronoiTesselation::dimension |
( |
| ) |
const |
|
inline |
Gets the dimension of the points.
Can be smaller than the dimension of the mesh.
Definition at line 308 of file CVT.h.
◆ done_current()
void GEO::CentroidalVoronoiTesselation::done_current |
( |
| ) |
|
|
inline |
◆ embedding()
double* GEO::CentroidalVoronoiTesselation::embedding |
( |
index_t |
p | ) |
|
|
inline |
Returns the representation of a point in embedding space.
- Parameters
-
- Returns
- a pointer to the coordinates of the point
- Precondition
- p < nb_points()
Definition at line 335 of file CVT.h.
◆ funcgrad()
virtual void GEO::CentroidalVoronoiTesselation::funcgrad |
( |
index_t |
n, |
|
|
double * |
x, |
|
|
double & |
f, |
|
|
double * |
g |
|
) |
| |
|
protectedvirtual |
Computes the objective function and its gradient.
- Parameters
-
[in] | n | number of variables |
[in] | x | current value of the variables |
[out] | f | current value of the objective function |
[out] | g | gradient of the objective function |
◆ funcgrad_CB()
static void GEO::CentroidalVoronoiTesselation::funcgrad_CB |
( |
index_t |
n, |
|
|
double * |
x, |
|
|
double & |
f, |
|
|
double * |
g |
|
) |
| |
|
static |
Callback for the numerical solver.
Evaluates the objective function and its gradient.
- Parameters
-
[in] | n | number of variables |
[in] | x | current value of the variables |
[out] | f | current value of the objective function |
[out] | g | gradient of the objective function |
◆ Lloyd_iterations()
virtual void GEO::CentroidalVoronoiTesselation::Lloyd_iterations |
( |
index_t |
nb_iter | ) |
|
|
virtual |
Relaxes the points with Lloyd's algorithm.
It is in general less efficient than Newton, but more resistant to heterogeneous point distribution. Therefore a small number of Lloyd iterations may be used right after a call to compute_initial_sampling() to regularize the point set before calling Newton_iterations().
- Parameters
-
[in] | nb_iter | number of iterations |
◆ lock_point()
void GEO::CentroidalVoronoiTesselation::lock_point |
( |
index_t |
i | ) |
|
|
inline |
Locks a point.
A locked point is constrained to stay at the same position during the optimization.
- Parameters
-
- Precondition
- i < nb_points()
Definition at line 377 of file CVT.h.
◆ make_current()
void GEO::CentroidalVoronoiTesselation::make_current |
( |
| ) |
|
|
inline |
◆ mesh()
Mesh* GEO::CentroidalVoronoiTesselation::mesh |
( |
| ) |
|
|
inline |
Returns the input mesh.
Definition at line 213 of file CVT.h.
◆ newiteration()
virtual void GEO::CentroidalVoronoiTesselation::newiteration |
( |
| ) |
|
|
protectedvirtual |
Callback for the numerical solver.
Updates the progress bar.
◆ newiteration_CB()
static void GEO::CentroidalVoronoiTesselation::newiteration_CB |
( |
index_t |
n, |
|
|
const double * |
x, |
|
|
double |
f, |
|
|
const double * |
g, |
|
|
double |
gnorm |
|
) |
| |
|
static |
Callback for the numerical solver.
Updates the progress bar.
- Parameters
-
[in] | n | number of variables |
[in] | x | current value of the variables |
[in] | f | current value of the objective function |
[in] | g | gradient of the objective function |
[in] | gnorm | norm of the gradient of the objective function |
◆ Newton_iterations()
virtual void GEO::CentroidalVoronoiTesselation::Newton_iterations |
( |
index_t |
nb_iter, |
|
|
index_t |
m = 7 |
|
) |
| |
|
virtual |
Relaxes the points with Newton-Lloyd's algorithm.
- Parameters
-
[in] | nb_iter | number of iterations |
[in] | m | number of evaluations used for Hessian approximation |
◆ point_is_locked()
bool GEO::CentroidalVoronoiTesselation::point_is_locked |
( |
index_t |
i | ) |
const |
|
inline |
Tests whether a point is locked.
A locked point is constrained to stay at the same position during the optimization.
- Parameters
-
- Precondition
- i < nb_points()
Definition at line 363 of file CVT.h.
◆ R3_embedding()
const vec3& GEO::CentroidalVoronoiTesselation::R3_embedding |
( |
index_t |
p | ) |
const |
|
inline |
Gets the representation of a point in R3.
- Parameters
-
- Returns
- a const reference to the 3d version of the point
- Precondition
- p < nb_points()
Definition at line 325 of file CVT.h.
◆ resize_points()
void GEO::CentroidalVoronoiTesselation::resize_points |
( |
index_t |
nb_points | ) |
|
Changes the number of points.
- Parameters
-
[in] | nb_points | new number of points |
Resizes the internal vector used to store the points
◆ RVD()
◆ set_constrained_cvt()
void GEO::CentroidalVoronoiTesselation::set_constrained_cvt |
( |
bool |
x | ) |
|
|
inline |
Specifies whether constrained mode should be used.
- Parameters
-
[in] | x | If set (default = false), compute_surface() projects the vertices onto the input surface. |
Definition at line 206 of file CVT.h.
◆ set_facets_range()
void GEO::CentroidalVoronoiTesselation::set_facets_range |
( |
index_t |
facets_begin, |
|
|
index_t |
facets_end |
|
) |
| |
|
inline |
Restricts computation to a part of the input mesh.
The part of the input mesh should be specified as a contiguous range of facet indices.
- Parameters
-
[in] | facets_begin | first facet in the range |
[in] | facets_end | one past last facet in the range |
Definition at line 238 of file CVT.h.
◆ set_points()
void GEO::CentroidalVoronoiTesselation::set_points |
( |
index_t |
nb_points, |
|
|
const double * |
points |
|
) |
| |
Initializes the points with a user-specified vector.
- Parameters
-
[in] | nb_points | number of points in points |
[in] | points | (size = dimension()*nb_points): user-defined initialization. It is copied into an internal vector |
◆ set_progress_logger()
void GEO::CentroidalVoronoiTesselation::set_progress_logger |
( |
ProgressTask * |
progress | ) |
|
|
inline |
Sets a client for the progress bars.
- Parameters
-
Definition at line 300 of file CVT.h.
◆ set_show_iterations()
void GEO::CentroidalVoronoiTesselation::set_show_iterations |
( |
bool |
x | ) |
|
|
inline |
Specifies whether a progress bar should be used.
- Parameters
-
[in] | x | If set, shows iterations using a "progress bar". |
Definition at line 187 of file CVT.h.
◆ set_use_RVC_centroids()
void GEO::CentroidalVoronoiTesselation::set_use_RVC_centroids |
( |
bool |
x | ) |
|
|
inline |
Specifies whether centroids of Voronoi cells should be used.
- Parameters
-
[in] | x | If set (default = true), compute_surface() replaces the vertices with the centroids of the connected components of the restricted Voronoi cells. |
Definition at line 197 of file CVT.h.
◆ set_volumetric()
void GEO::CentroidalVoronoiTesselation::set_volumetric |
( |
bool |
x | ) |
|
|
inline |
Sets volumetric mode.
- Parameters
-
[in] | x | if true, volumetric mode is used, otherwise surfacic mode is used. |
Definition at line 352 of file CVT.h.
◆ unlock_all_points()
void GEO::CentroidalVoronoiTesselation::unlock_all_points |
( |
| ) |
|
|
inline |
Unlocks all the points.
A locked point is constrained to stay at the same position during the optimization.
Definition at line 407 of file CVT.h.
◆ unlock_point()
void GEO::CentroidalVoronoiTesselation::unlock_point |
( |
index_t |
i | ) |
|
|
inline |
Unlocks a point.
A locked point is constrained to stay at the same position during the optimization.
- Parameters
-
- Precondition
- i < nb_points()
Definition at line 392 of file CVT.h.
◆ is_projection_
bool GEO::CentroidalVoronoiTesselation::is_projection_ |
|
protected |
the Nd -> 3d transform is a projection
Definition at line 457 of file CVT.h.
The documentation for this class was generated from the following file: