CGAL 4.9 - Scale-Space Surface Reconstruction
|
#include <CGAL/Scale_space_surface_reconstruction_3.h>
computes a triangulated surface mesh interpolating a point set.
This class stores several of the (intermediate) results. This makes it easier and more efficient to adjust the parameter settings based on preliminary results, or to further increase the scale to improve the results. The class stores the point set at the current scale and the reconstructed surface, possibly with iterators over the shells.
The class also stores the parameters for estimating the optimal neighborhood radius and either the lastest estimate or the manually set radius. This way, the radius can be estimated (again) whenever necessary. Also note that both increasing the scale and reconstructing the surface use this radius. By changing or re-estimating the radius between these operations, they can use separate parameter settings.
The surface can be constructed either for a fixed neighborhood radius, or for a dynamic radius. When constructing the surface for exactly one neighborhood radius, it is faster to set FS
to Tag_true
. If the correct neighborhood radius should be changed or estimated multiple times, it is faster to set FS
to Tag_false
.
It is undefined whether a surface with fixed radius may have its radius changed, but if so, this will likely require more computation time than changing the radius of a dynamic surface. In either case, it is possible to change the point set while maintaining the same radius.
The surface can be stored either as an unordered collection of triangles, or as a collection ordered by shells. A shell is a maximally connected component of the surface where connected facets are locally oriented towards the same side of the surface.
Gt | is the geometric traits class. It must be a model of DelaunayTriangulationTraits_3 . It must have a RealEmbeddable field number type. Generally, Exact_predicates_inexact_constructions_kernel is preferred. |
FS | determines whether the surface is expected to be constructed for a fixed neighborhood radius. It must be a Boolean_tag type. The default value is Tag_true . Note that the value of this parameter does not change the result but only has an impact on the run-time. |
wA | must be a model of DiagonalizeTraits and determines how to diagonalize a weighted covariance matrix to approximate a weighted point set. It can be omitted: if Eigen 3 (or greater) is available and CGAL_EIGEN3_ENABLED is defined then an overload using Eigen_diagonalize_traits is provided. Otherwise, the internal implementation Diagonalize_traits is used. |
Ct | indicates whether to use concurrent processing. It must be either Sequential_tag or Parallel_tag (the default value). |
Public Types | |
typedef Gt::Point_3 | Point |
defines the point type. | |
Types | |
typedef Gt::FT | FT |
defines the field number type. | |
typedef Gt::Vector_3 | Vector |
defines the vector type. | |
typedef Gt::Plane_3 | Plane |
defines the plane type. | |
typedef Gt::Triangle_3 | Triangle |
defines the triangle type. | |
typedef unspecified_type | Point_iterator |
defines an iterator over the points. | |
typedef const unspecified_type | Point_const_iterator |
defines a constant iterator over the points. | |
typedef CGAL::cpp11::array < unsigned int, 3 > | Triple |
defines a triple of point indices indicating a triangle of the surface. | |
typedef unspecified_type | Triple_iterator |
defines an iterator over the triples. | |
typedef const unspecified_type | Triple_const_iterator |
defines a constant iterator over the triples. | |
Constructors | |
Scale_space_surface_reconstruction_3 (unsigned int neighbors, unsigned int samples) | |
constructs a surface reconstructor that will automatically estimate the neighborhood radius. More... | |
Scale_space_surface_reconstruction_3 (FT sq_radius) | |
constructs a surface reconstructor with a given neighborhood radius. More... | |
Point Set Manipulation | |
template<class InputIterator > | |
void | insert (InputIterator begin, InputIterator end) |
inserts a collection of points into the scale-space at the current scale. More... | |
void | insert (const Point &p) |
inserts a point into the scale-space at the current scale. More... | |
void | clear () |
clears the stored scale-space surface reconstruction data. More... | |
Neighborhood Size Estimation | |
FT | estimate_neighborhood_squared_radius () |
estimates the neighborhood radius. More... | |
FT | estimate_neighborhood_squared_radius (unsigned int neighbors, unsigned int samples) |
estimates the neighborhood radius based on a number of sample points. More... | |
void | set_neighborhood_squared_radius (const FT &sq_radius) |
sets the squared radius of the neighborhood. More... | |
FT | neighborhood_squared_radius () const |
gives the squared radius of the neighborhood. More... | |
bool | has_neighborhood_squared_radius () const |
checks whether the neighborhood radius has been set. More... | |
Neighborhood Size Estimation Parameters | |
unsigned int | mean_number_of_neighbors () const |
gives the mean number of neighbors an estimated neighborhood should contain. More... | |
unsigned int | neighborhood_sample_size () const |
gives the number of sample points the neighborhood estimation uses. More... | |
void | set_mean_number_of_neighbors (unsigned int neighbors) |
sets the mean number of neighbors an estimated neighborhood should contain. More... | |
void | set_neighborhood_sample_size (unsigned int samples) |
sets the number of sample points the neighborhood estimation uses. More... | |
Scale-Space Manipulation | |
void | increase_scale (unsigned int iterations=1) |
increases the scale by a number of iterations. More... | |
Surface Reconstruction | |
void | reconstruct_surface (bool separate_shells=true, bool force_manifold=false, FT border_angle=45) |
constructs a triangle mesh from the point set at a fixed scale. More... | |
std::size_t | number_of_triangles () const |
gives the number of triangles of the surface. | |
std::size_t | number_of_points () const |
gives the number of points of the surface. | |
std::size_t | number_of_shells () const |
gives the number of shells of the surface. | |
Iterators | |
Point_const_iterator | points_begin () const |
gives an iterator to the first point at the current scale. | |
Point_iterator | points_begin () |
gives an iterator to the first point at the current scale. More... | |
Point_const_iterator | points_end () const |
gives a past-the-end iterator of the points at the current scale. | |
Point_iterator | points_end () |
gives a past-the-end iterator of the points at the current scale. More... | |
Triple_const_iterator | surface_begin () const |
gives an iterator to the first triple in the surface. | |
Triple_iterator | surface_begin () |
gives an iterator to the first triple in the surface. More... | |
Triple_const_iterator | surface_end () const |
gives a past-the-end iterator of the triples in the surface. | |
Triple_iterator | surface_end () |
gives a past-the-end iterator of the triples in the surface. More... | |
Triple_const_iterator | shell_begin (std::size_t shell) const |
gives an iterator to the first triple in a given shell. More... | |
Triple_iterator | shell_begin (std::size_t shell) |
gives an iterator to the first triple in a given shell. More... | |
Triple_const_iterator | shell_end (std::size_t shell) const |
gives a past-the-end iterator of the triples in a given shell. More... | |
Triple_iterator | shell_end (std::size_t shell) |
gives a past-the-end iterator of the triples in a given shell. More... | |
Triple_const_iterator | garbage_begin () const |
gives an iterator to the first triple of the garbage facets that may be discarded if 2-manifold output is required. More... | |
Triple_iterator | garbage_begin () |
gives an iterator to the first triple of the garbage facets that may be discarded if 2-manifold output is required. More... | |
Triple_const_iterator | garbage_end () const |
gives a past-the-end iterator of the triples of the garbage facets that may be discarded if 2-manifold output is required. More... | |
Triple_iterator | garbage_end () |
gives a past-the-end iterator of the triples of the garbage facets that may be discarded if 2-manifold output is required. More... | |
CGAL::Scale_space_surface_reconstruction_3< Gt, FS, wA, Ct >::Scale_space_surface_reconstruction_3 | ( | unsigned int | neighbors, |
unsigned int | samples | ||
) |
constructs a surface reconstructor that will automatically estimate the neighborhood radius.
neighbors | is the number of neighbors a point's neighborhood should contain on average. |
samples | is the number of points sampled to estimate the neighborhood radius. |
CGAL::Scale_space_surface_reconstruction_3< Gt, FS, wA, Ct >::Scale_space_surface_reconstruction_3 | ( | FT | sq_radius) |
constructs a surface reconstructor with a given neighborhood radius.
sq_radius | is the squared radius of the neighborhood. |
sq_radius
is not negative. void CGAL::Scale_space_surface_reconstruction_3< Gt, FS, wA, Ct >::clear | ( | ) |
clears the stored scale-space surface reconstruction data.
This includes discarding the surface, the scale-space and all its points, and any estimation of the neighborhood radius.
Methods called after this point may have to re-estimate the neighborhood radius. This method does not discard the parameters for estimating this radius (the mean number of neighbors and the sample size).
FT CGAL::Scale_space_surface_reconstruction_3< Gt, FS, wA, Ct >::estimate_neighborhood_squared_radius | ( | ) |
estimates the neighborhood radius.
This method is equivalent to running estimate_neighborhood_squared_radius( mean_number_of_neighbors(), neighborhood_sample_size() )
.
This method can be called by the scale-space and surface construction methods if the neighborhood radius is not set when they are called.
insert(begin, end)
.Scale_space_surface_reconstruction_3< Gt, FS, wA, Ct >::FT CGAL::Scale_space_surface_reconstruction_3< Gt, FS, wA, Ct >::estimate_neighborhood_squared_radius | ( | unsigned int | neighbors, |
unsigned int | samples | ||
) |
estimates the neighborhood radius based on a number of sample points.
The neighborhood radius is expressed as the radius of the smallest ball centered on a point such that the ball contains at least a specified number of points, not counting the point itself.
The neighborhood radius is set to the mean of these radii, taken over a number of point samples.
neighbors | is the number of neighbors a point's neighborhood should contain on average, not counting the point itself. |
samples | is the number of points sampled to estimate the neighborhood radius. |
insert(begin, end)
.Triple_const_iterator CGAL::Scale_space_surface_reconstruction_3< Gt, FS, wA, Ct >::garbage_begin | ( | ) | const |
gives an iterator to the first triple of the garbage facets that may be discarded if 2-manifold output is required.
Triple_iterator CGAL::Scale_space_surface_reconstruction_3< Gt, FS, wA, Ct >::garbage_begin | ( | ) |
gives an iterator to the first triple of the garbage facets that may be discarded if 2-manifold output is required.
Triple_const_iterator CGAL::Scale_space_surface_reconstruction_3< Gt, FS, wA, Ct >::garbage_end | ( | ) | const |
gives a past-the-end iterator of the triples of the garbage facets that may be discarded if 2-manifold output is required.
Triple_iterator CGAL::Scale_space_surface_reconstruction_3< Gt, FS, wA, Ct >::garbage_end | ( | ) |
gives a past-the-end iterator of the triples of the garbage facets that may be discarded if 2-manifold output is required.
bool CGAL::Scale_space_surface_reconstruction_3< Gt, FS, wA, Ct >::has_neighborhood_squared_radius | ( | ) | const |
checks whether the neighborhood radius has been set.
The radius can be set manually, or estimated automatically.
true
iff the radius has been either set manually or estimated.void CGAL::Scale_space_surface_reconstruction_3< Gt, FS, wA, Ct >::increase_scale | ( | unsigned int | iterations = 1 ) |
increases the scale by a number of iterations.
Each iteration the scale is increased, the points set at a higher scale is computed. At a higher scale, the points set is smoother.
If the neighborhood radius has not been set before, it is automatically estimated using estimate_neighborhood_squared_radius()
.
iterations | is the number of iterations to perform. If iterations is 0, nothing happens. |
insert(begin, end)
.void CGAL::Scale_space_surface_reconstruction_3< Gt, FS, wA, Ct >::insert | ( | InputIterator | begin, |
InputIterator | end | ||
) |
inserts a collection of points into the scale-space at the current scale.
InputIterator | is an iterator over the point collection. The value type of the iterator must be a Point . |
begin | is an iterator to the first point of the collection. |
end | is a past-the-end iterator for the point collection. |
reconstruct_surface()
after inserting the points.insert(const Point& p)
. void CGAL::Scale_space_surface_reconstruction_3< Gt, FS, wA, Ct >::insert | ( | const Point & | p) |
inserts a point into the scale-space at the current scale.
p | is the point to insert. |
reconstruct_surface()
.unsigned int CGAL::Scale_space_surface_reconstruction_3< Gt, FS, wA, Ct >::mean_number_of_neighbors | ( | ) | const |
gives the mean number of neighbors an estimated neighborhood should contain.
This number is only used if the neighborhood radius has not been set manually.
When the neighborhood radius is estimated, it should on average contain this many neighbors, not counting the neighborhood center.
unsigned int CGAL::Scale_space_surface_reconstruction_3< Gt, FS, wA, Ct >::neighborhood_sample_size | ( | ) | const |
gives the number of sample points the neighborhood estimation uses.
This number is only used if the neighborhood radius has not been set manually.
If the number of samples is larger than the point cloud, every point is used and the optimal neighborhood radius is computed exactly instead of estimated.
FT CGAL::Scale_space_surface_reconstruction_3< Gt, FS, wA, Ct >::neighborhood_squared_radius | ( | ) | const |
gives the squared radius of the neighborhood.
The neighborhood radius is used by increase_scale()
to compute the point set at the desired scale and by reconstruct_surface()
to construct a surface from the point set at the current scale.
Point_iterator CGAL::Scale_space_surface_reconstruction_3< Gt, FS, wA, Ct >::points_begin | ( | ) |
gives an iterator to the first point at the current scale.
Point_iterator CGAL::Scale_space_surface_reconstruction_3< Gt, FS, wA, Ct >::points_end | ( | ) |
gives a past-the-end iterator of the points at the current scale.
void CGAL::Scale_space_surface_reconstruction_3< Gt, FS, wA, Ct >::reconstruct_surface | ( | bool | separate_shells = true , |
bool | force_manifold = false , |
||
FT | border_angle = 45 |
||
) |
constructs a triangle mesh from the point set at a fixed scale.
The order of the points at the current scale is the same as the order at the original scale, meaning that the surface can interpolate the point set at the original scale by applying the indices of the surface triangles to the original point set.
After construction, the triangles of the surface can be iterated using surface_begin()
and surface_end()
.
If the neighborhood radius has not been set before, it is automatically estimated using estimate_neighborhood_squared_radius()
.
separate_shells | determines whether to collect the surface per shell. |
force_manifold | determines if the surface is forced to be 2-manifold. |
border_angle | sets the maximal angle between two facets such that the edge is seen as a border. |
If the output is forced to be 2-manifold, some almost flat volume bubbles are detected. To do so, border edges must be estimated.
An edge adjacent to 2 regular facets is considered as a border if it is also adjacent to a singular facet or if the angle between the two regular facets is lower than this parameter (set to 45° by default).
insert(begin, end)
. border_angle
is not used if force_manifold
is set to false.void CGAL::Scale_space_surface_reconstruction_3< Gt, FS, wA, Ct >::set_mean_number_of_neighbors | ( | unsigned int | neighbors) |
sets the mean number of neighbors an estimated neighborhood should contain.
This number is only used if the neighborhood radius has not been set manually.
When the neighborhood radius is estimated, it should on average contain this many neighbors, not counting the neighborhood center.
neighbors | is the number of neighbors a neighborhood ball centered on a point should contain on average when the radius is estimated, not counting the point itself. |
void CGAL::Scale_space_surface_reconstruction_3< Gt, FS, wA, Ct >::set_neighborhood_sample_size | ( | unsigned int | samples) |
sets the number of sample points the neighborhood estimation uses.
This number is only used if the neighborhood radius has not been set manually.
If the number of samples is larger than the point cloud, every point is used and the optimal neighborhood radius is computed exactly instead of estimated.
samples | is the number of points to sample for neighborhood estimation. |
void CGAL::Scale_space_surface_reconstruction_3< Gt, FS, wA, Ct >::set_neighborhood_squared_radius | ( | const FT & | sq_radius) |
sets the squared radius of the neighborhood.
The neighborhood radius is used by #[increase_scale()
to compute the point set at the desired scale and by reconstruct_surface()
to construct a surface from the point set at the current scale.
sq_radius | is the squared radius of the neighborhood. |
Scale_space_surface_reconstruction_3< Gt, FS, wA, Ct >::Triple_const_iterator CGAL::Scale_space_surface_reconstruction_3< Gt, FS, wA, Ct >::shell_begin | ( | std::size_t | shell) | const |
gives an iterator to the first triple in a given shell.
shell | is the index of the shell to access. |
shell
is in the range [ 0, number_of_shells()
). Scale_space_surface_reconstruction_3< Gt, FS, wA, Ct >::Triple_iterator CGAL::Scale_space_surface_reconstruction_3< Gt, FS, wA, Ct >::shell_begin | ( | std::size_t | shell) |
gives an iterator to the first triple in a given shell.
shell | is the index of the shell to access. |
shell
is in the range [ 0, number_of_shells()
).Scale_space_surface_reconstruction_3< Gt, FS, wA, Ct >::Triple_const_iterator CGAL::Scale_space_surface_reconstruction_3< Gt, FS, wA, Ct >::shell_end | ( | std::size_t | shell) | const |
gives a past-the-end iterator of the triples in a given shell.
shell | is the index of the shell to access. |
shell
is in the range [ 0, number_of_shells()
). Scale_space_surface_reconstruction_3< Gt, FS, wA, Ct >::Triple_iterator CGAL::Scale_space_surface_reconstruction_3< Gt, FS, wA, Ct >::shell_end | ( | std::size_t | shell) |
gives a past-the-end iterator of the triples in a given shell.
shell | is the index of the shell to access. |
shell
is in the range [ 0, number_of_shells()
).Triple_iterator CGAL::Scale_space_surface_reconstruction_3< Gt, FS, wA, Ct >::surface_begin | ( | ) |
gives an iterator to the first triple in the surface.
Triple_iterator CGAL::Scale_space_surface_reconstruction_3< Gt, FS, wA, Ct >::surface_end | ( | ) |
gives a past-the-end iterator of the triples in the surface.