DGtal 1.3.0
Loading...
Searching...
No Matches
Data Structures | Public Types | Public Member Functions | Static Public Attributes
DGtal::QuickHull< TKernel > Struct Template Reference

Aim: Implements the quickhull algorithm by Barber et al. [9], a famous arbitrary dimensional convex hull computation algorithm. It relies on dedicated geometric kernels for computing and comparing facet geometries. More...

#include <DGtal/geometry/tools/QuickHull.h>

Data Structures

struct  Facet
 

Public Types

enum  { UNASSIGNED = (Index) -1 }
 Label for points that are not assigned to any facet. More...
 
enum class  Status {
  Uninitialized = 0 , InputInitialized = 1 , SimplexCompleted = 2 , FacetsCompleted = 3 ,
  VerticesCompleted = 4 , AllCompleted = 5 , NotFullDimensional = 10 , InvalidRidge = 11 ,
  InvalidConvexHull = 12
}
 Represents the status of a QuickHull object. More...
 
typedef TKernel Kernel
 
typedef Kernel::CoordinatePoint Point
 
typedef Kernel::CoordinateVector Vector
 
typedef Kernel::CoordinateScalar Scalar
 
typedef Kernel::InternalScalar InternalScalar
 
typedef std::size_t Index
 
typedef std::size_t Size
 
typedef std::vector< IndexIndexRange
 
typedef Kernel::HalfSpace HalfSpace
 
typedef Kernel::CombinatorialPlaneSimplex CombinatorialPlaneSimplex
 
typedef std::pair< Index, IndexRidge
 

Public Member Functions

 BOOST_STATIC_ASSERT ((Point::dimension==Vector::dimension))
 
Standard services (construction, initialization, accessors)
 QuickHull (const Kernel &K=Kernel(), int dbg=0)
 
Status status () const
 
void clear ()
 Clears the object. More...
 
Size memory () const
 
Size nbPoints () const
 
Size nbFacets () const
 
Size nbVertices () const
 
Size nbFiniteFacets () const
 
Size nbInfiniteFacets () const
 
Initialization services
template<typename InputPoint >
bool setInput (const std::vector< InputPoint > &input_points, bool remove_duplicates=true)
 
bool setInitialSimplex (const IndexRange &full_splx)
 
Convex hull services
bool computeConvexHull (Status target=Status::VerticesCompleted)
 
bool computeInitialSimplex ()
 
bool computeFacets ()
 
bool computeVertices ()
 
Output services
template<typename OutputPoint >
bool getVertexPositions (std::vector< OutputPoint > &vertex_positions)
 
bool getVertex2Point (IndexRange &vertex_to_point)
 
bool getPoint2Vertex (IndexRange &point_to_vertex)
 
bool getFacetVertices (std::vector< IndexRange > &facet_vertices) const
 
bool getFacetHalfSpaces (std::vector< HalfSpace > &facet_halfspaces)
 
Check hull services
bool check ()
 
bool checkFacets ()
 
bool checkHull ()
 
Interface
void selfDisplay (std::ostream &out) const
 
bool isValid () const
 

Data Fields

public datas
Kernel kernel
 Kernel that is duplicated for computing facet geometries. More...
 
int debug_level
 debug_level from 0:no to 2 More...
 
std::vector< Pointpoints
 the set of points, indexed as in the array. More...
 
IndexRange input2comp
 
IndexRange comp2input
 
IndexRange processed_points
 Points already processed (and within the convex hull). More...
 
std::vector< Indexassignment
 assignment of points to facets More...
 
std::vector< Facetfacets
 the current set of facets. More...
 
std::set< Indexdeleted_facets
 set of deleted facets More...
 
IndexRange p2v
 point index -> vertex index (or UNASSIGNED) More...
 
IndexRange v2p
 vertex index -> point index More...
 
Size nb_finite_facets
 Number of finite facets. More...
 
Size nb_infinite_facets
 Number of infinite facets (!= 0 only for specific kernels) More...
 
std::vector< double > timings
 Timings of the different phases: 0: init, 1: facets, 2: vertices. More...
 
std::vector< Sizefacet_counter
 Counts the number of facets with a given number of vertices. More...
 

Static Public Attributes

static const Size dimension = Point::dimension
 

Protected Attributes

protected datas
Status myStatus
 

protected services

InternalScalar height (const Facet &F, const Point &p) const
 
bool above (const Facet &F, const Point &p) const
 
bool aboveOrOn (const Facet &F, const Point &p) const
 
bool on (const Facet &F, const Point &p) const
 
void cleanFacets ()
 
void cleanInfiniteFacets ()
 
void renumberInfiniteFacets ()
 
bool processFacet (std::queue< Index > &Q)
 
bool checkFacet (Index f) const
 
Index newFacet ()
 
void deleteFacet (Index f)
 
void makeNeighbors (const Index if1, const Index if2)
 
void unmakeNeighbors (const Index if1, const Index if2)
 
Index mergeParallelFacets (const Index if1, const Index if2)
 
bool areFacetsParallel (const Index if1, const Index if2) const
 
bool areFacetsNeighbor (const Index if1, const Index if2) const
 
Facet makeFacet (const IndexRange &base, Index below) const
 
IndexRange pointsOnRidge (const Ridge &R) const
 
IndexRange orientedFacetPoints (Index f) const
 
void filterVerticesOnFacet (const Index f)
 
IndexRange pickInitialSimplex () const
 
bool computeSimplexConfiguration (const IndexRange &full_simplex)
 
static IndexRange pickIntegers (const Size d, const Size n)
 

Detailed Description

template<typename TKernel>
struct DGtal::QuickHull< TKernel >

Aim: Implements the quickhull algorithm by Barber et al. [9], a famous arbitrary dimensional convex hull computation algorithm. It relies on dedicated geometric kernels for computing and comparing facet geometries.

Description of template class 'QuickHull'

You can use it to build convex hulls of points with integral coordinate (using kernel ConvexHullIntegralKernel), convex hulls with rational coordinates (using kernel ConvexHullRationalKernel) or to build Delaunay triangulations (using kernels DelaunayIntegralKernel and DelaunayRationalKernel).

See also
QuickHull algorithm in arbitrary dimension for convex hull and Delaunay cell complex computation

Below is a complete example that computes the convex hull of points randomly defined in a ball, builds a 3D mesh out of it and output it as an OBJ file.

#include "DGtal/base/Common.h"
#include "DGtal/kernel/PointVector.h"
#include "DGtal/shapes/SurfaceMesh.h"
#include "DGtal/io/writers/SurfaceMeshWriter.h"
#include "QuickHull.h"
using namespace DGtal::Z3i;
int main( int argc, char* argv[] )
{
int nb = argc > 1 ? atoi( argv[ 1 ] ) : 100; // nb points
int R = argc > 2 ? atoi( argv[ 2 ] ) : 10; // x-radius of ellipsoid
// (0) typedefs
// (1) create range of random points in ball
std::vector< Point > V;
const auto R2 = R
for ( int i = 0; i < nb; ) {
Point p( rand() % (2*R+1) - R, rand() % (2*R+1) - R, rand() % (2*R+1) - R );
if ( p.squaredNorm() < R2 ) { V.push_back( p ); i++; }
}
// (2) compute convex hull
hull.setInput( V );
std::cout << "#points=" << hull.nbPoints()
<< " #vertices=" << hull.nbVertices()
<< " #facets=" << hull.nbFacets() << std::endl;
// (3) build mesh
std::vector< RealPoint > positions;
hull.getVertexPositions( positions );
std::vector< std::vector< std::size_t > > facets;
SMesh mesh( positions.cbegin(), positions.cend(), facets.cbegin(), facets.cend() );
// (4) output result as OBJ file
std::ofstream out( "qhull.obj" );
::writeOBJ( out, mesh );
out.close();
return 0;
}
SurfaceMesh< RealPoint, RealVector > SMesh
Z3i this namespace gathers the standard of types for 3D imagery.
Aim: a geometric kernel to compute the convex hull of digital points with integer-only arithmetic.
Aim: Implements the quickhull algorithm by Barber et al. , a famous arbitrary dimensional convex hull...
Definition: QuickHull.h:140
bool setInput(const std::vector< InputPoint > &input_points, bool remove_duplicates=true)
Definition: QuickHull.h:382
std::vector< Facet > facets
the current set of facets.
Definition: QuickHull.h:814
bool getFacetVertices(std::vector< IndexRange > &facet_vertices) const
Definition: QuickHull.h:666
Kernel::CoordinatePoint Point
Definition: QuickHull.h:142
Size nbVertices() const
Definition: QuickHull.h:337
bool getVertexPositions(std::vector< OutputPoint > &vertex_positions)
Definition: QuickHull.h:612
Size nbFacets() const
Definition: QuickHull.h:328
bool computeConvexHull(Status target=Status::VerticesCompleted)
Definition: QuickHull.h:460
Size nbPoints() const
Definition: QuickHull.h:323
Aim: An helper class for writing mesh file formats (Waverfront OBJ at this point) and creating a Surf...
Aim: Represents an embedded mesh as faces and a list of vertices. Vertices may be shared among faces ...
Definition: SurfaceMesh.h:92
int main()
Definition: testBits.cpp:56
Note
In opposition with the usual QuickHull implementation, this class uses a kernel that can be chosen in order to provide exact computations. This is the case for lattice points.
In opposition with CGAL 3D convex hull package, or with the arbitrary dimensional dD Triangulation package, this algorithm does not build a simplicial convex hull. Facets may not be trangles or simplices in higher dimensions.
This version is generally more than twice faster than CGAL convex_hull_3 for the usual CGAL kernels Cartesian and Exact_predicates_inexact_constructions_kernel.
However this implementation is not tailored for incremental dynamic convex hull computations.
Template Parameters
TKernelany type of QuickHull kernel, like ConvexHullIntegralKernel.
Examples
geometry/tools/checkLatticeBallQuickHull.cpp, geometry/tools/exampleLatticeBallDelaunay2D.cpp, geometry/tools/exampleLatticeBallQuickHull3D.cpp, and geometry/tools/exampleRationalBallQuickHull3D.cpp.

Definition at line 139 of file QuickHull.h.

Member Typedef Documentation

◆ CombinatorialPlaneSimplex

template<typename TKernel >
typedef Kernel::CombinatorialPlaneSimplex DGtal::QuickHull< TKernel >::CombinatorialPlaneSimplex

Definition at line 151 of file QuickHull.h.

◆ HalfSpace

template<typename TKernel >
typedef Kernel::HalfSpace DGtal::QuickHull< TKernel >::HalfSpace

Definition at line 150 of file QuickHull.h.

◆ Index

template<typename TKernel >
typedef std::size_t DGtal::QuickHull< TKernel >::Index

Definition at line 146 of file QuickHull.h.

◆ IndexRange

template<typename TKernel >
typedef std::vector< Index > DGtal::QuickHull< TKernel >::IndexRange

Definition at line 149 of file QuickHull.h.

◆ InternalScalar

template<typename TKernel >
typedef Kernel::InternalScalar DGtal::QuickHull< TKernel >::InternalScalar

Definition at line 145 of file QuickHull.h.

◆ Kernel

template<typename TKernel >
typedef TKernel DGtal::QuickHull< TKernel >::Kernel

Definition at line 141 of file QuickHull.h.

◆ Point

template<typename TKernel >
typedef Kernel::CoordinatePoint DGtal::QuickHull< TKernel >::Point

Definition at line 142 of file QuickHull.h.

◆ Ridge

template<typename TKernel >
typedef std::pair< Index, Index > DGtal::QuickHull< TKernel >::Ridge

A ridge for point p is a pair of facets, such that p is visible from first facet, but not from second facet.

Definition at line 236 of file QuickHull.h.

◆ Scalar

template<typename TKernel >
typedef Kernel::CoordinateScalar DGtal::QuickHull< TKernel >::Scalar

Definition at line 144 of file QuickHull.h.

◆ Size

template<typename TKernel >
typedef std::size_t DGtal::QuickHull< TKernel >::Size

Definition at line 147 of file QuickHull.h.

◆ Vector

template<typename TKernel >
typedef Kernel::CoordinateVector DGtal::QuickHull< TKernel >::Vector

Definition at line 143 of file QuickHull.h.

Member Enumeration Documentation

◆ anonymous enum

template<typename TKernel >
anonymous enum

Label for points that are not assigned to any facet.

Enumerator
UNASSIGNED 

Definition at line 155 of file QuickHull.h.

155{ UNASSIGNED = (Index) -1 };
std::size_t Index
Definition: QuickHull.h:146

◆ Status

template<typename TKernel >
enum class DGtal::QuickHull::Status
strong

Represents the status of a QuickHull object.

Enumerator
Uninitialized 

QuickHull is empty and has just been instantiated.

InputInitialized 

A range of input points has been given to QuickHull.

SimplexCompleted 

An initial full-dimensional simplex has been found. QuickHull core algorithm can start.

FacetsCompleted 

All facets of the convex hull are identified.

VerticesCompleted 

All vertices of the convex hull are determined.

AllCompleted 

Same as VerticesCompleted.

NotFullDimensional 

Error: the initial simplex is not full dimensional.

InvalidRidge 

Error: some ridge is not consistent (probably integer overflow).

InvalidConvexHull 

Error: the convex hull does not contain all its points (probably integer overflow).

Definition at line 239 of file QuickHull.h.

239 {
240 Uninitialized = 0,
241 InputInitialized = 1,
242 SimplexCompleted = 2,
243 FacetsCompleted = 3,
245 AllCompleted = 5,
247 InvalidRidge = 11,
249 };
@ InvalidRidge
Error: some ridge is not consistent (probably integer overflow).
@ InvalidConvexHull
Error: the convex hull does not contain all its points (probably integer overflow).
@ SimplexCompleted
An initial full-dimensional simplex has been found. QuickHull core algorithm can start.
@ FacetsCompleted
All facets of the convex hull are identified.
@ NotFullDimensional
Error: the initial simplex is not full dimensional.
@ VerticesCompleted
All vertices of the convex hull are determined.
@ InputInitialized
A range of input points has been given to QuickHull.
@ AllCompleted
Same as VerticesCompleted.
@ Uninitialized
QuickHull is empty and has just been instantiated.

Constructor & Destructor Documentation

◆ QuickHull()

template<typename TKernel >
DGtal::QuickHull< TKernel >::QuickHull ( const Kernel K = Kernel(),
int  dbg = 0 
)
inline

Default constructor

Parameters
[in]Ka kernel for computing facet geometries.
[in]dbgthe trace level, from 0 (no) to 3 (very verbose).

Definition at line 259 of file QuickHull.h.

261 {}
Kernel kernel
Kernel that is duplicated for computing facet geometries.
Definition: QuickHull.h:798
int debug_level
debug_level from 0:no to 2
Definition: QuickHull.h:800
KSpace K

Member Function Documentation

◆ above()

template<typename TKernel >
bool DGtal::QuickHull< TKernel >::above ( const Facet F,
const Point p 
) const
inlineprotected

◆ aboveOrOn()

template<typename TKernel >
bool DGtal::QuickHull< TKernel >::aboveOrOn ( const Facet F,
const Point p 
) const
inlineprotected
Parameters
Fany valid facet
pany point
Returns
'true' iff p is above F or on F.

Definition at line 862 of file QuickHull.h.

863 { return kernel.aboveOrOn( F.H, p ); }

References DGtal::QuickHull< TKernel >::Facet::H, and DGtal::QuickHull< TKernel >::kernel.

Referenced by DGtal::QuickHull< TKernel >::checkFacet(), and DGtal::QuickHull< TKernel >::processFacet().

◆ areFacetsNeighbor()

template<typename TKernel >
bool DGtal::QuickHull< TKernel >::areFacetsNeighbor ( const Index  if1,
const Index  if2 
) const
inlineprotected

Checks if two facets are neighbors by looking at the points on their boundary.

Parameters
if1a valid facet index
if2a valid facet index
Returns
'true' if the two facets have enough common points to be direct neighbors.

Definition at line 1292 of file QuickHull.h.

1293 {
1294 const Facet& f1 = facets[ if1 ];
1295 const Facet& f2 = facets[ if2 ];
1296 Index nb = 0;
1297 for ( Index i1 = 0, i2 = 0; i1 < f1.on_set.size() && i2 < f2.on_set.size(); )
1298 {
1299 if ( f1.on_set[ i1 ] == f2.on_set[ i2 ] ) { nb++; i1++; i2++; }
1300 else if ( f1.on_set[ i1 ] < f2.on_set[ i2 ] ) i1++;
1301 else i2++;
1302 }
1303 return nb >= ( dimension - 1 );
1304 }
SMesh::Index Index
static const Size dimension
Definition: QuickHull.h:152

References DGtal::QuickHull< TKernel >::dimension, DGtal::QuickHull< TKernel >::facets, and DGtal::QuickHull< TKernel >::Facet::on_set.

Referenced by DGtal::QuickHull< TKernel >::checkFacet(), and DGtal::QuickHull< TKernel >::processFacet().

◆ areFacetsParallel()

template<typename TKernel >
bool DGtal::QuickHull< TKernel >::areFacetsParallel ( const Index  if1,
const Index  if2 
) const
inlineprotected

Checks if two distinct facets are parallel (i.e. should be merged).

Parameters
if1a valid facet index
if2a valid facet index
Returns
'true' if the two facets are parallel.

Definition at line 1273 of file QuickHull.h.

1274 {
1275 ASSERT( if1 != if2 );
1276 const Facet& f1 = facets[ if1 ];
1277 const Facet& f2 = facets[ if2 ];
1278 if ( kernel.equal( f1.H, f2.H ) ) return true;
1279 // Need to check if one N is a multiple of the other.
1280 for ( auto&& v : f1.on_set )
1281 if ( ! on( f2, points[ v ] ) ) return false;
1282 return true;
1283 }
bool on(const Facet &F, const Point &p) const
Definition: QuickHull.h:868

References DGtal::QuickHull< TKernel >::facets, DGtal::QuickHull< TKernel >::Facet::H, DGtal::QuickHull< TKernel >::kernel, DGtal::QuickHull< TKernel >::on(), and DGtal::QuickHull< TKernel >::Facet::on_set.

Referenced by DGtal::QuickHull< TKernel >::processFacet().

◆ BOOST_STATIC_ASSERT()

template<typename TKernel >
DGtal::QuickHull< TKernel >::BOOST_STATIC_ASSERT ( (Point::dimension==Vector::dimension)  )

◆ check()

template<typename TKernel >
bool DGtal::QuickHull< TKernel >::check ( )
inline

Global validity check of the convex hull after processing.

Note
Be careful, this function is much slower than computing the convex hull. It can take a long time since its complexity is \( O(n f ) \), where n is the number of input points and f the number of facets.

Definition at line 715 of file QuickHull.h.

716 {
717 bool ok = true;
719 trace.warning() << "[Quickhull::check] invalid status="
720 << (int)status() << std::endl;
721 ok = false;
722 }
723 if ( processed_points.size() != points.size() ) {
724 trace.warning() << "[Quickhull::check] not all points processed: "
725 << processed_points.size() << " / " << points.size()
726 << std::endl;
727 ok = false;
728 }
729 if ( ! checkHull() ) {
730 trace.warning() << "[Quickhull::check] Check hull is invalid. "
731 << std::endl;
732 ok = false;
733 }
734 if ( ! checkFacets() ) {
735 trace.warning() << "[Quickhull::check] Check facets is invalid. "
736 << std::endl;
737 ok = false;
738 }
739 return ok;
740 }
std::ostream & warning()
Trace trace
Definition: Common.h:154
IndexRange processed_points
Points already processed (and within the convex hull).
Definition: QuickHull.h:810
bool checkFacets()
Definition: QuickHull.h:743
Status status() const
Definition: QuickHull.h:266
bool checkHull()
Definition: QuickHull.h:762

References DGtal::QuickHull< TKernel >::AllCompleted, DGtal::QuickHull< TKernel >::checkFacets(), DGtal::QuickHull< TKernel >::checkHull(), DGtal::QuickHull< TKernel >::FacetsCompleted, DGtal::QuickHull< TKernel >::processed_points, DGtal::QuickHull< TKernel >::status(), DGtal::trace, and DGtal::Trace::warning().

◆ checkFacet()

template<typename TKernel >
bool DGtal::QuickHull< TKernel >::checkFacet ( Index  f) const
inlineprotected
Returns
true if the facet is valid

Definition at line 1161 of file QuickHull.h.

1162 {
1163 const Facet& F = facets[ f ];
1164 bool ok = F.on_set.size() >= dimension;
1165 for ( auto v : F.on_set )
1166 if ( ! on( F, points[ v ] ) ) {
1167 trace.error() << "[QuickHull::checkFacet( " << f
1168 << ") Invalid 'on' vertex " << v << std::endl;
1169 ok = false;
1170 }
1171 if ( F.neighbors.size() < dimension ) {
1172 trace.error() << "[QuickHull::checkFacet( " << f
1173 << ") Not enough neighbors " << F.neighbors.size() << std::endl;
1174 ok = false;
1175 }
1176 for ( auto nf : F.neighbors )
1177 if ( ! areFacetsNeighbor( f, nf ) ) {
1178 trace.error() << "[QuickHull::checkFacet( " << f
1179 << ") Invalid neighbor " << nf << std::endl;
1180 ok = false;
1181 }
1182 if ( aboveOrOn( F, points[ F.below ] ) ) {
1183 trace.error() << "[QuickHull::checkFacet( " << f
1184 << ") Bad below point " << F.below << std::endl;
1185 ok = false;
1186 }
1187 for ( auto ov : F.outside_set )
1188 if ( ! above( F, points[ ov ] ) ) {
1189 trace.error() << "[QuickHull::checkFacet( " << f
1190 << ") Bad outside point " << ov << std::endl;
1191 ok = false;
1192 }
1193 for ( auto && v : F.on_set ) {
1194 Size n = 0;
1195 for ( auto&& N : facets[ f ].neighbors )
1196 if ( on( facets[ N ], points[ v ] ) ) n += 1;
1197 if ( n < dimension-1 ) {
1198 trace.error() << "[QuickHull::checkFacet( " << f << ") 'on' point " << v
1199 << " is a vertex of " << n << " facets "
1200 << "(should be >= " << dimension-1 << ")" << std::endl;
1201 ok = false;
1202 }
1203 }
1204 return ok;
1205 }
std::ostream & error()
bool areFacetsNeighbor(const Index if1, const Index if2) const
Definition: QuickHull.h:1292
bool aboveOrOn(const Facet &F, const Point &p) const
Definition: QuickHull.h:862
bool above(const Facet &F, const Point &p) const
Definition: QuickHull.h:856
HalfEdgeDataStructure::Size Size

References DGtal::QuickHull< TKernel >::above(), DGtal::QuickHull< TKernel >::aboveOrOn(), DGtal::QuickHull< TKernel >::areFacetsNeighbor(), DGtal::QuickHull< TKernel >::Facet::below, DGtal::QuickHull< TKernel >::dimension, DGtal::Trace::error(), DGtal::QuickHull< TKernel >::facets, DGtal::QuickHull< TKernel >::Facet::neighbors, DGtal::QuickHull< TKernel >::on(), DGtal::QuickHull< TKernel >::Facet::on_set, DGtal::QuickHull< TKernel >::Facet::outside_set, and DGtal::trace.

Referenced by DGtal::QuickHull< TKernel >::checkFacets(), and DGtal::QuickHull< TKernel >::processFacet().

◆ checkFacets()

template<typename TKernel >
bool DGtal::QuickHull< TKernel >::checkFacets ( )
inline
Returns
'true' if all facets are consistent with their neighors

Definition at line 743 of file QuickHull.h.

744 {
745 Size nb = 0;
746 Size nbok = 0;
747 for ( Index f = 0; f < facets.size(); ++f )
748 if ( deleted_facets.count( f ) == 0 ) {
749 bool ok = checkFacet( f );
750 nbok += ok ? 1 : 0;
751 nb += 1;
752 }
753 return nb == nbok;
754 }
bool checkFacet(Index f) const
Definition: QuickHull.h:1161
std::set< Index > deleted_facets
set of deleted facets
Definition: QuickHull.h:816

References DGtal::QuickHull< TKernel >::checkFacet(), DGtal::QuickHull< TKernel >::deleted_facets, and DGtal::QuickHull< TKernel >::facets.

Referenced by DGtal::QuickHull< TKernel >::check(), DGtal::QuickHull< TKernel >::computeSimplexConfiguration(), and DGtal::QuickHull< TKernel >::processFacet().

◆ checkHull()

template<typename TKernel >
bool DGtal::QuickHull< TKernel >::checkHull ( )
inline
Returns
'true' iff the current set of facets encloses all the points
Note
Be careful, this function is much slower than computing the convex hull. It can take a long time since its complexity is \( O(n f ) \), where n is the number of input points and f the number of facets.

Definition at line 762 of file QuickHull.h.

763 {
764 Index nbok = 0;
765 // for ( Index v = 0; v < points.size(); ++v ) {
766 for ( auto v : processed_points ) {
767 bool ok = true;
768 for ( Index f = 0; f < facets.size(); ++f )
769 if ( deleted_facets.count( f ) == 0 ) {
770 if ( above( facets[ f ], points[ v ] ) ) {
771 ok = false;
772 trace.error() << "- bad vertex " << v << " " << points[ v ]
773 << " dist="
774 << ( kernel.height( facets[ f ].H, points[ v ] ) )
775 << std::endl;
776 break;
777 }
778 }
779 nbok += ok ? 1 : 0;
780 }
781 if ( debug_level >= 2 ) {
782 trace.info() << nbok << "/"
783 << processed_points.size() << " vertices inside convex hull."
784 << std::endl;
785 }
787 return nbok == processed_points.size();
788 }
std::ostream & info()

References DGtal::QuickHull< TKernel >::above(), DGtal::QuickHull< TKernel >::debug_level, DGtal::QuickHull< TKernel >::deleted_facets, DGtal::Trace::error(), DGtal::QuickHull< TKernel >::facets, DGtal::H, DGtal::Trace::info(), DGtal::QuickHull< TKernel >::InvalidConvexHull, DGtal::QuickHull< TKernel >::kernel, DGtal::QuickHull< TKernel >::myStatus, DGtal::QuickHull< TKernel >::processed_points, and DGtal::trace.

Referenced by DGtal::QuickHull< TKernel >::check(), DGtal::QuickHull< TKernel >::computeSimplexConfiguration(), and DGtal::QuickHull< TKernel >::processFacet().

◆ cleanFacets()

template<typename TKernel >
void DGtal::QuickHull< TKernel >::cleanFacets ( )
inlineprotected

Cleans and renumber the facets so that no one belongs to deleted_facets.

Definition at line 873 of file QuickHull.h.

874 {
875 if ( deleted_facets.empty() ) return;
876 IndexRange renumbering( facets.size() );
877 Index i = 0;
878 Index j = 0;
879 for ( auto& l : renumbering ) {
880 if ( ! deleted_facets.count( j ) ) l = i++;
881 else l = UNASSIGNED;
882 j++;
883 }
884 const Index nf = facets.size() - deleted_facets.size();
885 deleted_facets.clear();
886 for ( Index f = 0; f < facets.size(); f++ )
887 if ( ( renumbering[ f ] != UNASSIGNED ) && ( f != renumbering[ f ] ) )
888 facets[ renumbering[ f ] ] = facets[ f ];
889 facets.resize( nf );
890 for ( auto& F : facets ) {
891 for ( auto& N : F.neighbors ) {
892 if ( renumbering[ N ] == UNASSIGNED )
893 trace.error() << "Invalid deleted neighboring facet." << std::endl;
894 else N = renumbering[ N ];
895 }
896 }
897 }
std::vector< Index > IndexRange
Definition: QuickHull.h:149

References DGtal::QuickHull< TKernel >::deleted_facets, DGtal::Trace::error(), DGtal::QuickHull< TKernel >::facets, DGtal::trace, and DGtal::QuickHull< TKernel >::UNASSIGNED.

Referenced by DGtal::QuickHull< TKernel >::computeFacets().

◆ cleanInfiniteFacets()

template<typename TKernel >
void DGtal::QuickHull< TKernel >::cleanInfiniteFacets ( )
inlineprotected

Removes infinite facets and renumber the facets in a consecutive way.

Deprecated:

Definition at line 901 of file QuickHull.h.

902 {
903 if ( ! kernel.hasInfiniteFacets() ) return;
904 IndexRange renumbering( facets.size() );
905 Index i = 0;
906 Index j = 0;
907 for ( auto& l : renumbering ) {
908 if ( ! kernel.isHalfSpaceFacetInfinite( facets[ j ].H ) ) l = i++;
909 else l = UNASSIGNED;
910 j++;
911 }
912 const Index nf = i;
913 for ( Index f = 0; f < facets.size(); f++ )
914 if ( ( renumbering[ f ] != UNASSIGNED ) && ( f != renumbering[ f ] ) )
915 facets[ renumbering[ f ] ] = facets[ f ];
916 facets.resize( nf );
917 for ( auto& F : facets ) {
918 for ( auto& N : F.neighbors ) {
919 N = renumbering[ N ];
920 }
921 }
922 }

References DGtal::QuickHull< TKernel >::facets, DGtal::QuickHull< TKernel >::kernel, and DGtal::QuickHull< TKernel >::UNASSIGNED.

◆ clear()

template<typename TKernel >
void DGtal::QuickHull< TKernel >::clear ( )
inline

Clears the object.

Definition at line 270 of file QuickHull.h.

271 {
273 points.clear();
274 processed_points.clear();
275 input2comp.clear();
276 comp2input.clear();
277 assignment.clear();
278 facets.clear();
279 deleted_facets.clear();
280 p2v.clear();
281 v2p.clear();
282 timings.clear();
283 }
IndexRange v2p
vertex index -> point index
Definition: QuickHull.h:820
IndexRange input2comp
Definition: QuickHull.h:805
std::vector< Index > assignment
assignment of points to facets
Definition: QuickHull.h:812
IndexRange p2v
point index -> vertex index (or UNASSIGNED)
Definition: QuickHull.h:818
IndexRange comp2input
Definition: QuickHull.h:808
std::vector< double > timings
Timings of the different phases: 0: init, 1: facets, 2: vertices.
Definition: QuickHull.h:826

References DGtal::QuickHull< TKernel >::assignment, DGtal::QuickHull< TKernel >::comp2input, DGtal::QuickHull< TKernel >::deleted_facets, DGtal::QuickHull< TKernel >::facets, DGtal::QuickHull< TKernel >::input2comp, DGtal::QuickHull< TKernel >::myStatus, DGtal::QuickHull< TKernel >::p2v, DGtal::QuickHull< TKernel >::processed_points, DGtal::QuickHull< TKernel >::timings, DGtal::QuickHull< TKernel >::Uninitialized, and DGtal::QuickHull< TKernel >::v2p.

Referenced by DGtal::QuickHull< TKernel >::setInput().

◆ computeConvexHull()

template<typename TKernel >
bool DGtal::QuickHull< TKernel >::computeConvexHull ( Status  target = Status::VerticesCompleted)
inline

Computes the convex hull of the given range of points until a specified target:

  • Status::SimplexCompleted: an initial full dimensional simplex is determined.
  • Status::FacetsCompleted: all the facets of the convex hull are extracted (core of quickhull). You can stop here if you only need an H-representation of the output polytope.
  • Status::VerticesCompleted: all the vertices of the hull are determined and are consistently ordered on each facet. You need to stop here if you need a V-representation of the output polytope.
Precondition
status() must be at least Status::InputInitialized
Parameters
[in]targetthe computation target in Status::SimplexCompleted, Status::FacetsCompleted, Status::VerticesCompleted.
Returns
'true' if the computation target has been successfully achieved, 'false' if the achieved status is not the one specified.
Examples
geometry/tools/exampleLatticeBallQuickHull3D.cpp.

Definition at line 460 of file QuickHull.h.

461 {
462 if ( target < Status::InputInitialized || target > Status::AllCompleted )
463 return false;
464 Clock tic;
466 { // Initialization
467 tic.startClock();
468 bool ok1 = computeInitialSimplex();
469 timings.push_back( tic.stopClock() );
470 if ( ! ok1 ) return false;
471 if ( status() == target ) return true;
472 }
474 { // Computes facets
475 tic.startClock();
476 bool ok2 = computeFacets();
477 timings.push_back( tic.stopClock() );
478 if ( ! ok2 ) return false;
479 if ( status() == target ) return true;
480 }
482 { // Computes vertices
483 tic.startClock();
484 bool ok3 = computeVertices();
485 timings.push_back( tic.stopClock() );
486 if ( ! ok3 ) return false;
487 if ( status() == target ) return true;
488 }
489 if ( target == Status::AllCompleted
491 { // for now, Status::VerticesCompleted and
492 // Status::AllCompleted are the same.
494 return true;
495 }
496 return false;
497 }
void tic()
Starts timer.
bool computeInitialSimplex()
Definition: QuickHull.h:504
bool computeFacets()
Definition: QuickHull.h:523
bool computeVertices()
Definition: QuickHull.h:553

References DGtal::QuickHull< TKernel >::AllCompleted, DGtal::QuickHull< TKernel >::computeFacets(), DGtal::QuickHull< TKernel >::computeInitialSimplex(), DGtal::QuickHull< TKernel >::computeVertices(), DGtal::QuickHull< TKernel >::FacetsCompleted, DGtal::QuickHull< TKernel >::InputInitialized, DGtal::QuickHull< TKernel >::myStatus, DGtal::QuickHull< TKernel >::SimplexCompleted, DGtal::QuickHull< TKernel >::status(), tic(), DGtal::QuickHull< TKernel >::timings, and DGtal::QuickHull< TKernel >::VerticesCompleted.

Referenced by main().

◆ computeFacets()

template<typename TKernel >
bool DGtal::QuickHull< TKernel >::computeFacets ( )
inline

Computes the facets of the convex hull using Quickhull algorithm. If everything went well, the status is Status::FacetsCompleted afterwards.

Precondition
the status shoud be Status::SimplexCompleted (computeInitialSimplex should have been called).
Returns
'true' except if the status is not Initialized when called.

Definition at line 523 of file QuickHull.h.

524 {
525 if ( status() != Status::SimplexCompleted ) return false;
526 std::queue< Index > Q;
527 for ( Index fi = 0; fi < facets.size(); ++fi )
528 Q.push( fi );
529 Index n = 0;
530 while ( processFacet( Q ) ) {
531 if ( debug_level >= 1 )
532 trace.info() << "---- Iteration " << n++ << " #Q=" << Q.size() << std::endl;
533 }
534 cleanFacets();
535 if ( debug_level >= 2 ) {
536 trace.info() << ".... #facets=" << facets.size()
537 << " #deleted=" << deleted_facets.size() << std::endl;
538 }
540 return true;
541 }
bool processFacet(std::queue< Index > &Q)
Definition: QuickHull.h:969
void cleanFacets()
Definition: QuickHull.h:873

References DGtal::QuickHull< TKernel >::cleanFacets(), DGtal::QuickHull< TKernel >::debug_level, DGtal::QuickHull< TKernel >::deleted_facets, DGtal::QuickHull< TKernel >::facets, DGtal::QuickHull< TKernel >::FacetsCompleted, DGtal::Trace::info(), DGtal::QuickHull< TKernel >::myStatus, DGtal::QuickHull< TKernel >::processFacet(), DGtal::QuickHull< TKernel >::SimplexCompleted, DGtal::QuickHull< TKernel >::status(), and DGtal::trace.

Referenced by DGtal::QuickHull< TKernel >::computeConvexHull().

◆ computeInitialSimplex()

template<typename TKernel >
bool DGtal::QuickHull< TKernel >::computeInitialSimplex ( )
inline

Computes the initial full dimensional simplex from the input data.

Returns
'true' iff the input data contains d+1 points in general position, the object has then the status SimplexCompleted, otherwise returns 'false' and the status is NotFullDimensional.

Definition at line 504 of file QuickHull.h.

505 {
506 const auto full_simplex = pickInitialSimplex();
507 if ( full_simplex.empty() ) {
509 return false;
510 }
511 return computeSimplexConfiguration( full_simplex );
512 }
bool computeSimplexConfiguration(const IndexRange &full_simplex)
Definition: QuickHull.h:1454
IndexRange pickInitialSimplex() const
Definition: QuickHull.h:1397

References DGtal::QuickHull< TKernel >::computeSimplexConfiguration(), DGtal::QuickHull< TKernel >::myStatus, DGtal::QuickHull< TKernel >::NotFullDimensional, and DGtal::QuickHull< TKernel >::pickInitialSimplex().

Referenced by DGtal::QuickHull< TKernel >::computeConvexHull().

◆ computeSimplexConfiguration()

template<typename TKernel >
bool DGtal::QuickHull< TKernel >::computeSimplexConfiguration ( const IndexRange full_simplex)
inlineprotected

Computes the initial configuration induced by the given full dimensional simplex. Follows Status::InputInitialized and and terminate with Status::SimplexCompleted if everything went well.

Returns
'true' iff the input data contains d+1 points in general position, the object has then the status SimplexCompleted, otherwise returns 'false' and the status is NotFullDimensional.

Definition at line 1454 of file QuickHull.h.

1455 {
1456 assignment = std::vector< Index >( points.size(), UNASSIGNED );
1457 facets.resize( dimension + 1 );
1458 deleted_facets.clear();
1459 for ( Index j = 0; j < full_simplex.size(); ++j )
1460 {
1461 IndexRange lsimplex( dimension );
1462 IndexRange isimplex( dimension );
1463 Index s = 0;
1464 for ( Index i = 0; i <= dimension; i++ )
1465 if ( i != j ) {
1466 lsimplex[ s ] = i;
1467 isimplex[ s ] = full_simplex[ i ];
1468 s++;
1469 }
1470 facets[ j ] = makeFacet( isimplex, full_simplex[ j ] );
1471 facets[ j ].neighbors = lsimplex;
1472 for ( auto&& v : isimplex ) facets[ j ].on_set.push_back( v );
1473 std::sort( facets[ j ].on_set.begin(), facets[ j ].on_set.end() );
1474 }
1475 // List of unassigned vertices
1476 for ( Index fi = 0; fi < facets.size(); ++fi ) {
1477 Facet& f = facets[ fi ];
1478 for ( Index v = 0; v < points.size(); v++ )
1479 {
1480 if ( assignment[ v ] == UNASSIGNED && above( f, points[ v ] ) ) {
1481 f.outside_set.push_back( v );
1482 assignment[ v ] = fi;
1483 }
1484 }
1485 }
1486 for ( Index v = 0; v < points.size(); v++ )
1487 if ( assignment[ v ] == UNASSIGNED )
1488 processed_points.push_back( v );
1489
1490 // Display some information
1491 if ( debug_level >= 2 ) {
1492 for ( auto&& f : facets ) f.display( trace.info() );
1493 }
1495 if ( debug_level >= 1 ) {
1496 trace.info() << "[CHECK INVARIANT] " << processed_points.size()
1497 << " / " << points.size() << " points processed." << std::endl;
1498 bool okh = checkHull();
1499 if ( ! okh )
1500 trace.error() << "[computeInitialSimplex] Invalid convex hull" << std::endl;
1501 bool okf = checkFacets();
1502 if ( ! okf )
1503 trace.error() << "[computeInitialSimplex] Invalid facets" << std::endl;
1504 if ( ! ( okh && okf ) ) myStatus = Status::InvalidConvexHull;
1505 }
1506 return status() == Status::SimplexCompleted;
1507 }
Facet makeFacet(const IndexRange &base, Index below) const
Definition: QuickHull.h:1317

References DGtal::QuickHull< TKernel >::above(), DGtal::QuickHull< TKernel >::assignment, DGtal::QuickHull< TKernel >::checkFacets(), DGtal::QuickHull< TKernel >::checkHull(), DGtal::QuickHull< TKernel >::debug_level, DGtal::QuickHull< TKernel >::deleted_facets, DGtal::QuickHull< TKernel >::dimension, DGtal::Trace::error(), DGtal::QuickHull< TKernel >::facets, DGtal::Trace::info(), DGtal::QuickHull< TKernel >::InvalidConvexHull, DGtal::QuickHull< TKernel >::makeFacet(), DGtal::QuickHull< TKernel >::myStatus, DGtal::QuickHull< TKernel >::processed_points, DGtal::QuickHull< TKernel >::SimplexCompleted, DGtal::QuickHull< TKernel >::status(), DGtal::trace, and DGtal::QuickHull< TKernel >::UNASSIGNED.

Referenced by DGtal::QuickHull< TKernel >::computeInitialSimplex(), and DGtal::QuickHull< TKernel >::setInitialSimplex().

◆ computeVertices()

template<typename TKernel >
bool DGtal::QuickHull< TKernel >::computeVertices ( )
inline

Computes the vertices of the convex hull once the facets have been computed. It computes for each facet its vertices and reorder them so that, taken in order, their orientation matches the orientation of the facet. If everything went well, the status is Status::VerticesCompleted afterwards.

Precondition
the status shoud be Status::FacetsCompleted (computeFacets should have been called).
Returns
'true' except if the status is not FacetsCompleted when called.

Definition at line 553 of file QuickHull.h.

554 {
555 static const int MAX_NB_VPF = 10 * dimension;
556 if ( status() != Status::FacetsCompleted ) return false;
557
558 // Renumber infinite facets in case of Delaunay triangulation computation.
560
561 // Builds the maps v2p: vertex -> point, and p2v : point -> vertex.
562 facet_counter = IndexRange( MAX_NB_VPF, 0 );
563 v2p.clear();
564 p2v.resize( points.size() );
565 std::vector< IndexRange > p2f( points.size() );
566 for ( Index f = 0; f < facets.size(); ++f ) {
567 for ( auto&& p : facets[ f ].on_set ) p2f[ p ].push_back( f );
568 }
569
570 // vertices belong to at least d facets
571 Index v = 0;
572 for ( Index p = 0; p < points.size(); ++p ) {
573 const auto nbf = p2f[ p ].size();
574 facet_counter[ std::min( (int) nbf, MAX_NB_VPF-1 ) ] += 1;
575 if ( nbf >= dimension ) {
576 v2p.push_back( p );
577 p2v[ p ] = v++;
578 }
579 else p2v[ p ] = UNASSIGNED;
580 }
581
582 // Display debug informations
583 if ( debug_level >= 1 ) {
584 trace.info() << "#vertices=" << v2p.size() << " #facets=" << facets.size()
585 << std::endl;
586 trace.info() << "#inc_facets/point= ";
587 for ( auto n : facet_counter ) trace.info() << n << " ";
588 trace.info() << std::endl;
589 }
591 return true;
592 }
std::vector< Size > facet_counter
Counts the number of facets with a given number of vertices.
Definition: QuickHull.h:828
void renumberInfiniteFacets()
Definition: QuickHull.h:926

References DGtal::QuickHull< TKernel >::debug_level, DGtal::QuickHull< TKernel >::dimension, DGtal::QuickHull< TKernel >::facet_counter, DGtal::QuickHull< TKernel >::facets, DGtal::QuickHull< TKernel >::FacetsCompleted, DGtal::Trace::info(), DGtal::QuickHull< TKernel >::myStatus, DGtal::QuickHull< TKernel >::p2v, DGtal::QuickHull< TKernel >::renumberInfiniteFacets(), DGtal::QuickHull< TKernel >::status(), DGtal::trace, DGtal::QuickHull< TKernel >::UNASSIGNED, DGtal::QuickHull< TKernel >::v2p, and DGtal::QuickHull< TKernel >::VerticesCompleted.

Referenced by DGtal::QuickHull< TKernel >::computeConvexHull().

◆ deleteFacet()

template<typename TKernel >
void DGtal::QuickHull< TKernel >::deleteFacet ( Index  f)
inlineprotected

Deletes the given facet f.

Parameters
fa valid facet index

Definition at line 1218 of file QuickHull.h.

1219 {
1220 for ( auto n : facets[ f ].neighbors )
1221 facets[ n ].subNeighbor( f );
1222 deleted_facets.insert( f );
1223 facets[ f ].clear();
1224 }

References DGtal::QuickHull< TKernel >::deleted_facets, and DGtal::QuickHull< TKernel >::facets.

Referenced by DGtal::QuickHull< TKernel >::processFacet().

◆ filterVerticesOnFacet()

template<typename TKernel >
void DGtal::QuickHull< TKernel >::filterVerticesOnFacet ( const Index  f)
inlineprotected

Filters each vertex on the facet f to keep only the ones that are on or below the neighboring facets

Note
intended for debugging purposes.
Parameters
[in]fany valid facet index.

Definition at line 1377 of file QuickHull.h.

1378 {
1379 auto & on_set = facets[ f ].on_set;
1380 for ( Index i = 0; i < on_set.size(); )
1381 {
1382 Index v = on_set[ i ];
1383 Size n = 0;
1384 for ( auto&& N : facets[ f ].neighbors )
1385 if ( on( facets[ N ], points[ v ] ) ) n += 1;
1386 if ( n >= dimension-1 ) i++;
1387 else {
1388 on_set[ i ] = on_set.back();
1389 on_set.pop_back();
1390 }
1391 }
1392 std::sort( on_set.begin(), on_set.end() );
1393 }

References DGtal::QuickHull< TKernel >::dimension, DGtal::QuickHull< TKernel >::facets, and DGtal::QuickHull< TKernel >::on().

◆ getFacetHalfSpaces()

template<typename TKernel >
bool DGtal::QuickHull< TKernel >::getFacetHalfSpaces ( std::vector< HalfSpace > &  facet_halfspaces)
inline

This methods return the halfspaces corresponding to each facet. It is useful to build a polytope.

Parameters
[out]facet_halfspacesthe range giving for each facet its halfspace as a pair (normal, intercept).
Returns
'true' if the convex hull was computed before and status() was Status::FacetsCompleted, Status::VerticesCompleted or Status::AllCompleted, 'false' otherwise.

Definition at line 690 of file QuickHull.h.

691 {
692 facet_halfspaces.clear();
693 if ( ! ( status() >= Status::FacetsCompleted
694 && status() <= Status::AllCompleted ) ) return false;
695 facet_halfspaces.reserve( nbFacets() );
696 for ( Index f = 0; f < nbFacets(); ++f ) {
697 facet_halfspaces.push_back( facets[ f ].H );
698 }
699 return true;
700 }

References DGtal::QuickHull< TKernel >::AllCompleted, DGtal::QuickHull< TKernel >::facets, DGtal::QuickHull< TKernel >::FacetsCompleted, DGtal::H, DGtal::QuickHull< TKernel >::nbFacets(), and DGtal::QuickHull< TKernel >::status().

◆ getFacetVertices()

template<typename TKernel >
bool DGtal::QuickHull< TKernel >::getFacetVertices ( std::vector< IndexRange > &  facet_vertices) const
inline
Parameters
[out]facet_verticesthe range giving for each facet the indices of its vertices.
Returns
'true' if the convex hull was computed before and status() was Status::VerticesCompleted or Status::AllCompleted, 'false' otherwise.
Examples
geometry/tools/exampleLatticeBallQuickHull3D.cpp.

Definition at line 666 of file QuickHull.h.

667 {
668 facet_vertices.clear();
670 && status() <= Status::AllCompleted ) ) return false;
671 facet_vertices.reserve( nbFacets() );
672 for ( Index f = 0; f < nbFacets(); ++f ) {
673 IndexRange ofacet = orientedFacetPoints( f );
674 for ( auto& v : ofacet ) v = p2v[ v ];
675 facet_vertices.push_back( ofacet );
676 }
677 return true;
678 }
IndexRange orientedFacetPoints(Index f) const
Definition: QuickHull.h:1348

References DGtal::QuickHull< TKernel >::AllCompleted, DGtal::QuickHull< TKernel >::nbFacets(), DGtal::QuickHull< TKernel >::orientedFacetPoints(), DGtal::QuickHull< TKernel >::p2v, DGtal::QuickHull< TKernel >::status(), and DGtal::QuickHull< TKernel >::VerticesCompleted.

Referenced by main().

◆ getPoint2Vertex()

template<typename TKernel >
bool DGtal::QuickHull< TKernel >::getPoint2Vertex ( IndexRange point_to_vertex)
inline

Gets for each point its index in the output range of vertices (or UNASSIGNED if the point is not part of the convex hull)

Parameters
[out]point_to_vertexthe range giving for each point its index in the output range of vertices (or UNASSIGNED if the point is not part of the convex hull)
Returns
'true' if the convex hull was computed before and status() was Status::VerticesCompleted or Status::AllCompleted, 'false' otherwise.

Definition at line 651 of file QuickHull.h.

652 {
653 point_to_vertex.clear();
655 && status() <= Status::AllCompleted ) ) return false;
656 point_to_vertex = p2v;
657 return true;
658 }

References DGtal::QuickHull< TKernel >::AllCompleted, DGtal::QuickHull< TKernel >::p2v, DGtal::QuickHull< TKernel >::status(), and DGtal::QuickHull< TKernel >::VerticesCompleted.

◆ getVertex2Point()

template<typename TKernel >
bool DGtal::QuickHull< TKernel >::getVertex2Point ( IndexRange vertex_to_point)
inline

Gets for each vertex its index in the input range of points.

Parameters
[out]vertex_to_pointthe range giving for each vertex its index in the input range of points.
Returns
'true' if the convex hull was computed before and status() was Status::VerticesCompleted or Status::AllCompleted, 'false' otherwise.

Definition at line 632 of file QuickHull.h.

633 {
634 vertex_to_point.clear();
636 && status() <= Status::AllCompleted ) ) return false;
637 vertex_to_point = v2p;
638 return true;
639 }

References DGtal::QuickHull< TKernel >::AllCompleted, DGtal::QuickHull< TKernel >::status(), DGtal::QuickHull< TKernel >::v2p, and DGtal::QuickHull< TKernel >::VerticesCompleted.

◆ getVertexPositions()

template<typename TKernel >
template<typename OutputPoint >
bool DGtal::QuickHull< TKernel >::getVertexPositions ( std::vector< OutputPoint > &  vertex_positions)
inline

Gets the positions of the convex hull vertices.

Template Parameters
OutputPointa model of point such that the Point datatype is convertible to it.
Parameters
[out]vertex_positionsthe range of vertex positions.
Returns
'true' if the convex hull was computed before and status() was Status::VerticesCompleted or Status::AllCompleted, 'false' otherwise.
Examples
geometry/tools/exampleLatticeBallQuickHull3D.cpp.

Definition at line 612 of file QuickHull.h.

613 {
614 vertex_positions.clear();
616 && status() <= Status::AllCompleted ) ) return false;
617 vertex_positions.resize( v2p.size() );
618 for ( Index i = 0; i < v2p.size(); i++ ) {
619 kernel.convertPointTo( points[ v2p[ i ] ], vertex_positions[ i ] );
620 }
621 return true;
622 }

References DGtal::QuickHull< TKernel >::AllCompleted, DGtal::QuickHull< TKernel >::kernel, DGtal::QuickHull< TKernel >::status(), DGtal::QuickHull< TKernel >::v2p, and DGtal::QuickHull< TKernel >::VerticesCompleted.

Referenced by main().

◆ height()

template<typename TKernel >
InternalScalar DGtal::QuickHull< TKernel >::height ( const Facet F,
const Point p 
) const
inlineprotected
Parameters
Fany valid facet
pany point
Returns
the height of p wrt F (0: on, >0: above ).

Definition at line 850 of file QuickHull.h.

851 { return kernel.height( F.H, p ); }

References DGtal::QuickHull< TKernel >::Facet::H, and DGtal::QuickHull< TKernel >::kernel.

◆ isValid()

template<typename TKernel >
bool DGtal::QuickHull< TKernel >::isValid ( ) const
inline

Checks the validity/consistency of the object.

Returns
'true' if the object is valid, 'false' otherwise.

Definition at line 1540 of file QuickHull.h.

1541 {
1542 return status() >= Status::Uninitialized
1544 }

References DGtal::QuickHull< TKernel >::AllCompleted, DGtal::QuickHull< TKernel >::status(), and DGtal::QuickHull< TKernel >::Uninitialized.

◆ makeFacet()

template<typename TKernel >
Facet DGtal::QuickHull< TKernel >::makeFacet ( const IndexRange base,
Index  below 
) const
inlineprotected

Builds a facet from a base convex set of at least d different points and a point below.

Parameters
[in]basea range containing at least d distinct points in general position.
[in]belowa point below the hyperplane containing simplex.
Returns
a facet object of normal and c parameter such as the points of simplex lies on the facet hyperplane while point below lies under.

Definition at line 1317 of file QuickHull.h.

1318 {
1320 for ( Size i = 0; i < dimension; i++ ) simplex[ i ] = base[ i ];
1321 auto plane = kernel.compute( points, simplex, below );
1322 return Facet( plane, below );
1323 }
Kernel::CombinatorialPlaneSimplex CombinatorialPlaneSimplex
Definition: QuickHull.h:151

References DGtal::QuickHull< TKernel >::dimension, and DGtal::QuickHull< TKernel >::kernel.

Referenced by DGtal::QuickHull< TKernel >::computeSimplexConfiguration(), and DGtal::QuickHull< TKernel >::processFacet().

◆ makeNeighbors()

template<typename TKernel >
void DGtal::QuickHull< TKernel >::makeNeighbors ( const Index  if1,
const Index  if2 
)
inlineprotected

Makes two distinct facets if1 and if2 as neighbors

Parameters
if1a valid facet index
if2a valid facet index

Definition at line 1229 of file QuickHull.h.

1230 {
1231 facets[ if1 ].addNeighbor( if2 );
1232 facets[ if2 ].addNeighbor( if1 );
1233 }

References DGtal::QuickHull< TKernel >::facets.

Referenced by DGtal::QuickHull< TKernel >::mergeParallelFacets(), and DGtal::QuickHull< TKernel >::processFacet().

◆ memory()

template<typename TKernel >
Size DGtal::QuickHull< TKernel >::memory ( ) const
inline
Returns
an estimation of the current memory occupation of this object.

Definition at line 287 of file QuickHull.h.

288 {
289 // int debug_level;
290 Size M = sizeof( kernel ) + sizeof( int );
291 // std::vector< Point > points;
292 M += sizeof( std::vector< Point > )
293 + points.capacity() * sizeof( Point );
294 M += sizeof( std::vector< Index > )
295 + processed_points.capacity() * sizeof( Index );
296 M += sizeof( std::vector< Index > )
297 + input2comp.capacity() * sizeof( Index );
298 M += sizeof( std::vector< Index > )
299 + comp2input.capacity() * sizeof( Index );
300 // std::vector< Index > assignment;
301 M += sizeof( std::vector< Index > )
302 + assignment.capacity() * sizeof( Index );
303 // std::vector< Facet > facets;
304 M += sizeof( std::vector< Facet > )
305 + facets.capacity() * sizeof( Facet );
306 for ( const auto& f : facets ) M += f.variableMemory();
307 // std::set< Index > deleted_facets;
308 M += sizeof( std::set< Index > )
309 + deleted_facets.size() * ( sizeof( Index ) + 2*sizeof(Index*) );
310 // IndexRange p2v;
311 M += sizeof( std::vector< Index > )
312 + p2v.capacity() * sizeof( Index );
313 // IndexRange v2p;
314 M += sizeof( std::vector< Index > )
315 + v2p.capacity() * sizeof( Index );
316 // std::vector< double > timings;
317 M += sizeof( std::vector< double > )
318 + timings.capacity() * sizeof( double );
319 return M;
320 }

References DGtal::QuickHull< TKernel >::assignment, DGtal::QuickHull< TKernel >::comp2input, DGtal::QuickHull< TKernel >::deleted_facets, DGtal::QuickHull< TKernel >::facets, DGtal::QuickHull< TKernel >::input2comp, DGtal::QuickHull< TKernel >::kernel, DGtal::QuickHull< TKernel >::p2v, DGtal::QuickHull< TKernel >::processed_points, DGtal::QuickHull< TKernel >::timings, and DGtal::QuickHull< TKernel >::v2p.

◆ mergeParallelFacets()

template<typename TKernel >
Index DGtal::QuickHull< TKernel >::mergeParallelFacets ( const Index  if1,
const Index  if2 
)
inlineprotected

Merge the two facets and return the index of the second one, which is the one deleted.

Parameters
if1a valid facet index
if2a valid facet index
Precondition
the two facets should be distinct and parallel

Definition at line 1250 of file QuickHull.h.

1251 {
1252 Facet& f1 = facets[ if1 ];
1253 Facet& f2 = facets[ if2 ];
1254 std::copy( f2.outside_set.cbegin(), f2.outside_set.cend(),
1255 std::back_inserter( f1.outside_set ) );
1256 IndexRange merge_idx;
1257 std::set_union( f1.on_set.cbegin(), f1.on_set.cend(),
1258 f2.on_set.cbegin(), f2.on_set.cend(),
1259 std::back_inserter( merge_idx ) );
1260 f1.on_set.swap( merge_idx );
1261 for ( auto && nf2 : f2.neighbors ) {
1262 if ( nf2 == if1 ) continue;
1263 facets[ nf2 ].subNeighbor( if2 );
1264 makeNeighbors( if1, nf2 );
1265 }
1266 return if2;
1267 }
void makeNeighbors(const Index if1, const Index if2)
Definition: QuickHull.h:1229

References DGtal::QuickHull< TKernel >::facets, DGtal::QuickHull< TKernel >::makeNeighbors(), DGtal::QuickHull< TKernel >::Facet::neighbors, DGtal::QuickHull< TKernel >::Facet::on_set, and DGtal::QuickHull< TKernel >::Facet::outside_set.

Referenced by DGtal::QuickHull< TKernel >::processFacet().

◆ nbFacets()

template<typename TKernel >
Size DGtal::QuickHull< TKernel >::nbFacets ( ) const
inline

◆ nbFiniteFacets()

template<typename TKernel >
Size DGtal::QuickHull< TKernel >::nbFiniteFacets ( ) const
inline
Returns
the number of finite facets of the convex hull.
Precondition
status() >= Status::VerticesCompleted

Definition at line 346 of file QuickHull.h.

347 {
350 return nb_finite_facets;
351 }
Size nb_finite_facets
Number of finite facets.
Definition: QuickHull.h:822

References DGtal::QuickHull< TKernel >::AllCompleted, DGtal::QuickHull< TKernel >::nb_finite_facets, DGtal::QuickHull< TKernel >::status(), and DGtal::QuickHull< TKernel >::VerticesCompleted.

◆ nbInfiniteFacets()

template<typename TKernel >
Size DGtal::QuickHull< TKernel >::nbInfiniteFacets ( ) const
inline
Returns
the number of infinite facets of the convex hull.
Precondition
status() >= Status::VerticesCompleted

Definition at line 355 of file QuickHull.h.

356 {
359 return nb_infinite_facets;
360 }
Size nb_infinite_facets
Number of infinite facets (!= 0 only for specific kernels)
Definition: QuickHull.h:824

References DGtal::QuickHull< TKernel >::AllCompleted, DGtal::QuickHull< TKernel >::nb_infinite_facets, DGtal::QuickHull< TKernel >::status(), and DGtal::QuickHull< TKernel >::VerticesCompleted.

◆ nbPoints()

template<typename TKernel >
Size DGtal::QuickHull< TKernel >::nbPoints ( ) const
inline
Returns
the number of points used as input.
Examples
geometry/tools/exampleLatticeBallQuickHull3D.cpp.

Definition at line 323 of file QuickHull.h.

324 { return points.size(); }

Referenced by main(), and DGtal::QuickHull< TKernel >::selfDisplay().

◆ nbVertices()

template<typename TKernel >
Size DGtal::QuickHull< TKernel >::nbVertices ( ) const
inline

◆ newFacet()

template<typename TKernel >
Index DGtal::QuickHull< TKernel >::newFacet ( )
inlineprotected
Returns
an unused facet index

Definition at line 1208 of file QuickHull.h.

1209 {
1210 // SLightly faster to postpone deletion of intermediate facets.
1211 const Index f = facets.size();
1212 facets.push_back( Facet() );
1213 return f;
1214 }

References DGtal::QuickHull< TKernel >::facets.

Referenced by DGtal::QuickHull< TKernel >::processFacet().

◆ on()

template<typename TKernel >
bool DGtal::QuickHull< TKernel >::on ( const Facet F,
const Point p 
) const
inlineprotected
Parameters
Fany valid facet
pany point
Returns
'true' iff p lies on F.

Definition at line 868 of file QuickHull.h.

869 { return kernel.on( F.H, p ); }

References DGtal::QuickHull< TKernel >::Facet::H, and DGtal::QuickHull< TKernel >::kernel.

Referenced by DGtal::QuickHull< TKernel >::areFacetsParallel(), DGtal::QuickHull< TKernel >::checkFacet(), and DGtal::QuickHull< TKernel >::filterVerticesOnFacet().

◆ orientedFacetPoints()

template<typename TKernel >
IndexRange DGtal::QuickHull< TKernel >::orientedFacetPoints ( Index  f) const
inlineprotected

Given a facet index f, return its points oriented consistently with respect to the normal.

Parameters
fany valid facet index
Returns
the range of point indices that support this facet, in a consistent ordering.
Note
the order of points is consistent if, picking any d-simplex in order in this range, their associated half-spaces have all the same orientation.

Definition at line 1348 of file QuickHull.h.

1349 {
1350 const Facet& F = facets[ f ];
1351 IndexRange result = F.on_set;
1352 // Sort a facet such that its points, taken in order, have
1353 // always the same orientation of the facet. More precisely,
1354 // facets span a `dimension-1` vector space. There are thus
1355 // dimension-2 fixed points, and the last ones (at least two)
1356 // may be reordered.
1358 for ( Dimension k = 0; k < dimension-2; k++ )
1359 splx[ k ] = result[ k ];
1360 std::sort( result.begin()+dimension-2, result.end(),
1361 [&] ( Index i, Index j )
1362 {
1363 splx[ dimension-2 ] = i;
1364 splx[ dimension-1 ] = j;
1365 const auto H = kernel.compute( points, splx );
1366 const auto orient = kernel.dot( F.H, H );
1367 return orient > 0;
1368 } );
1369 return result;
1370 }
DGtal::uint32_t Dimension
Definition: Common.h:137

References DGtal::QuickHull< TKernel >::dimension, DGtal::QuickHull< TKernel >::facets, and DGtal::QuickHull< TKernel >::Facet::on_set.

Referenced by DGtal::QuickHull< TKernel >::getFacetVertices().

◆ pickInitialSimplex()

template<typename TKernel >
IndexRange DGtal::QuickHull< TKernel >::pickInitialSimplex ( ) const
inlineprotected
Returns
a full dimensional simplex as a vector of d + 1 distinct indices of input points, or an empty vector if none was found.

Definition at line 1397 of file QuickHull.h.

1398 {
1399 const Size nb = points.size();
1400 if ( nb < dimension + 1 ) return IndexRange();
1401 IndexRange best = pickIntegers( dimension + 1, nb );
1403 for ( Index j = 0; j < dimension; ++j ) splx[ j ] = best[ j ];
1404 const auto first_H = kernel.compute( points, splx, best.back() );
1405 auto best_volume = kernel.volume ( first_H, points[ best.back() ] );
1406 const Size nbtries = std::min( (Size) 10, 1 + nb / 10 );
1407 const Size max_nbtries = std::max( (Size) 10, 2 * nb );
1408 for ( Size i = 0; i < max_nbtries; i++ )
1409 {
1410 IndexRange tmp = pickIntegers( dimension + 1, nb );
1411 for ( Index j = 0; j < dimension; ++j ) splx[ j ] = tmp[ j ];
1412 const auto tmp_H = kernel.compute( points, splx, tmp.back() );
1413 const auto tmp_volume = kernel.volume ( tmp_H, points[ tmp.back() ] );
1414 if ( best_volume < tmp_volume ) {
1415 if ( debug_level >= 1 ) {
1416 trace.info() << "(" << i << ")"
1417 << " new_volume = " << tmp_volume
1418 << " > " << best_volume << std::endl;
1419 }
1420 best = tmp;
1421 best_volume = tmp_volume;
1422 }
1423 if ( i >= nbtries && best_volume > 0 ) return best;
1424 }
1425 return IndexRange();
1426 }
static IndexRange pickIntegers(const Size d, const Size n)
Definition: QuickHull.h:1431

References DGtal::QuickHull< TKernel >::debug_level, DGtal::QuickHull< TKernel >::dimension, DGtal::Trace::info(), DGtal::QuickHull< TKernel >::kernel, DGtal::QuickHull< TKernel >::pickIntegers(), and DGtal::trace.

Referenced by DGtal::QuickHull< TKernel >::computeInitialSimplex().

◆ pickIntegers()

template<typename TKernel >
static IndexRange DGtal::QuickHull< TKernel >::pickIntegers ( const Size  d,
const Size  n 
)
inlinestaticprotected
Returns
a vector of d distinct integers in {0, 1, ..., n-1} randomly chosen.
Parameters
[in]dthe number of returned integers
[in]nthe range of possible integers {0, 1, ..., n-1}

Definition at line 1431 of file QuickHull.h.

1432 {
1433 IndexRange result( d );
1434 bool distinct = false;
1435 while ( ! distinct )
1436 {
1437 distinct = true;
1438 for ( Index i = 0; i < d; i++ ) result[ i ] = rand() % n;
1439 std::sort( result.begin(), result.end() );
1440 for ( Index i = 1; distinct && i < d; i++ )
1441 distinct = result[ i-1 ] != result[ i ];
1442 }
1443 return result;
1444 }

Referenced by DGtal::QuickHull< TKernel >::pickInitialSimplex().

◆ pointsOnRidge()

template<typename TKernel >
IndexRange DGtal::QuickHull< TKernel >::pointsOnRidge ( const Ridge R) const
inlineprotected
Parameters
[in]Ra ridge between two facets (a pair of facets).
Returns
the points lying on a ridge.

Definition at line 1327 of file QuickHull.h.

1328 {
1329 // Slightly faster to look at points instead at true geometry.
1330 IndexRange result;
1331 const Facet& f1 = facets[ R.first ];
1332 const Facet& f2 = facets[ R.second ];
1333 std::set_intersection( f1.on_set.cbegin(), f1.on_set.cend(),
1334 f2.on_set.cbegin(), f2.on_set.cend(),
1335 std::back_inserter( result ) );
1336 return result;
1337 }

References DGtal::QuickHull< TKernel >::facets, DGtal::QuickHull< TKernel >::Facet::on_set, and DGtal::R.

Referenced by DGtal::QuickHull< TKernel >::processFacet().

◆ processFacet()

template<typename TKernel >
bool DGtal::QuickHull< TKernel >::processFacet ( std::queue< Index > &  Q)
inlineprotected

Process top facet in queue Q as in Quickhull algorithm, and updates the queue accordingly (top facet poped, new facets pushed).

Parameters
[in,out]Qa queue of facet index to process. which is updated by the method.
Returns
'true' if there is still work to do, 'false' when finished
Note
Core of Quickhull algorithm.

Definition at line 969 of file QuickHull.h.

970 {
971 // If Q empty, we are done
972 if ( Q.empty() ) return false;
973 Index F = Q.front();
974 Q.pop();
975 // If F is already deleted, proceed to next in queue.
976 if ( deleted_facets.count( F ) ) return true;
977 // Take car of current facet.
978 const Facet& facet = facets[ F ];
979 if ( debug_level >= 3 ) {
980 trace.info() << "---------------------------------------------" << std::endl;
981 trace.info() << "---- ACTIVE FACETS---------------------------" << std::endl;
982 bool ok = true;
983 for ( Index i = 0; i < facets.size(); i++ )
984 if ( ! deleted_facets.count( i ) ) {
985 trace.info() << "- facet " << i << " ";
986 facets[ i ].display( trace.info() );
987 ok = ok && checkFacet( i );
988 }
989 if ( ! ok ) { // should not happen.
991 return false;
992 }
993 }
994 if ( debug_level >= 2 ) {
995 trace.info() << "---------------------------------------------" << std::endl;
996 trace.info() << "Processing facet " << F << " ";
997 facet.display( trace.info() );
998 }
999 if ( facet.outside_set.empty() ) return true;
1000 // Selects furthest vertex
1001 Index furthest_v = facet.outside_set[ 0 ];
1002 auto furthest_h = height( facet, points[ furthest_v ] );
1003 for ( Index v = 1; v < facet.outside_set.size(); v++ ) {
1004 auto h = height( facet, points[ facet.outside_set[ v ] ] );
1005 if ( h > furthest_h ) {
1006 furthest_h = h;
1007 furthest_v = facet.outside_set[ v ];
1008 }
1009 }
1010 const Point& p = points[ furthest_v ];
1011 // Extracts Visible facets V and Horizon Ridges H
1012 std::vector< Index > V; // visible facets
1013 std::set< Index > M; // marked facets (are in E or were in E)
1014 std::queue< Index > E; // queue to extract visible facets
1015 std::vector< Ridge > H; // visible facets
1016 E.push ( F );
1017 M.insert( F );
1018 while ( ! E.empty() ) {
1019 Index G = E.front(); E.pop();
1020 V.push_back( G );
1021 for ( auto& N : facets[ G ].neighbors ) {
1022 if ( aboveOrOn( facets[ N ], p ) ) {
1023 if ( M.count( N ) ) continue;
1024 E.push( N );
1025 } else {
1026 H.push_back( { G, N } );
1027 }
1028 M.insert( N );
1029 }
1030 } // while ( ! E.empty() )
1031 if ( debug_level >= 1 ) {
1032 trace.info() << "#Visible=" << V.size() << " #Horizon=" << H.size()
1033 << " furthest_v=" << furthest_v << std::endl;
1034 }
1035 // Create new facets
1036 IndexRange new_facets;
1037 // For each ridge R in H
1038 for ( Index i = 0; i < H.size(); i++ )
1039 {
1040 // Create a new facet from R and p
1041 IndexRange ridge = pointsOnRidge( H[ i ] );
1042 if ( debug_level >= 3 ) {
1043 trace.info() << "Ridge (" << H[i].first << "," << H[i].second << ") = {";
1044 for ( auto&& r : ridge ) trace.info() << " " << r;
1045 trace.info() << " } furthest_v=" << furthest_v << std::endl;
1046 }
1047 IndexRange base( 1 + ridge.size() );
1048 Index j = 0;
1049 base[ j++ ] = furthest_v;
1050 for ( auto&& v : ridge ) base[ j++ ] = v;
1051 if ( j < dimension ) {
1052 trace.error() << "Bad ridge between " << std::endl
1053 << "- facet " << H[i].first << " ";
1054 facets[ H[i].first ].display( trace.error() );
1055 trace.error() << "- facet " << H[i].second << " ";
1056 facets[ H[i].second ].display( trace.error() );
1057 }
1058 Index nf = newFacet();
1059 new_facets.push_back( nf );
1060 facets[ nf ] = makeFacet( base, facets[ H[i].first ].below );
1061 facets[ nf ].on_set = IndexRange { base.cbegin(), base.cend() };
1062 std::sort( facets[ nf ].on_set.begin(), facets[ nf ].on_set.end() );
1063 makeNeighbors( nf, H[ i ].second );
1064 if ( debug_level >= 3 ) {
1065 trace.info() << "* New facet " << nf << " ";
1066 facets[ nf ].display( trace.info() );
1067 }
1068 // Checks that the facet is not parallel to another in the Horizon
1069 for ( Index k = 0; k < new_facets.size() - 1; k++ )
1070 if ( areFacetsParallel( new_facets[ k ], nf ) ) {
1071 if ( debug_level >= 1 ) {
1072 trace.info() << "Facets " << new_facets[ k ] << " and " << nf
1073 << " are parallel => merge." << std::endl;
1074 }
1075 mergeParallelFacets( new_facets[ k ], nf );
1076 new_facets.pop_back();
1077 deleteFacet( nf );
1078 if ( debug_level >= 3 ) {
1079 facets[ new_facets[ k ] ].display( trace.info() );
1080 }
1081 }
1082 }
1083 // For each new facet
1084 for ( Index i = 0; i < new_facets.size(); i++ )
1085 { // link the new facet to its neighbors
1086 for ( Index j = i + 1; j < new_facets.size(); j++ )
1087 {
1088 const Index nfi = new_facets[ i ];
1089 const Index nfj = new_facets[ j ];
1090 if ( areFacetsNeighbor( nfi, nfj ) )
1091 makeNeighbors( nfi, nfj );
1092 }
1093 }
1094 // Extracts all outside points from visible facets V
1095 IndexRange outside_pts;
1096 for ( auto&& vf : V ) {
1097 for ( auto&& v : facets[ vf ].outside_set ) {
1098 if ( v != furthest_v ) {
1099 outside_pts.push_back( v );
1100 assignment[ v ] = UNASSIGNED;
1101 }
1102 }
1103 }
1104 // For each new facet F'
1105 for ( Index i = 0; i < new_facets.size(); i++ ) {
1106 Facet& Fp = facets[ new_facets[ i ] ];
1107 Index max_j = outside_pts.size();
1108 for ( Index j = 0; j < max_j; ) {
1109 const Index v = outside_pts[ j ];
1110 if ( above( Fp, points[ v ] ) ) {
1111 Fp.outside_set.push_back( v );
1112 assignment[ v ] = new_facets[ i ];
1113 outside_pts[ j ] = outside_pts.back();
1114 outside_pts.pop_back();
1115 max_j--;
1116 } else j++;
1117 }
1118 if ( debug_level >= 3 ) {
1119 trace.info() << "- New facet " << new_facets[ i ] << " ";
1120 Fp.display( trace.info() );
1121 }
1122 }
1123 // Update processed points
1124 processed_points.push_back( furthest_v );
1125 for ( auto v : outside_pts ) processed_points.push_back( v );
1126
1127 // Delete the facets in V
1128 for ( auto&& v : V ) {
1129 if ( debug_level >= 2 ) {
1130 trace.info() << "Delete facet " << v << " ";
1131 facets[ v ].display( trace.info() );
1132 }
1133 deleteFacet( v );
1134 }
1135
1136 // Add new facets to queue
1137 for ( Index i = 0; i < new_facets.size(); i++ )
1138 Q.push( new_facets[ i ] );
1139 if ( debug_level >= 1 ) {
1140 trace.info() << "#facets=" << facets.size()
1141 << " #deleted=" << deleted_facets.size() << std::endl;
1142 }
1143
1144 // Checks that everything is ok.
1145 if ( debug_level >= 1 ) {
1146 trace.info() << "[CHECK INVARIANT] " << processed_points.size()
1147 << " / " << points.size() << " points processed." << std::endl;
1148 bool okh = checkHull();
1149 if ( ! okh )
1150 trace.error() << "[computeFacet] Invalid convex hull" << std::endl;
1151 bool okf = checkFacets();
1152 if ( ! okf )
1153 trace.error() << "[computeFacet] Invalid facets" << std::endl;
1154 if ( ! ( okh && okf ) ) myStatus = Status::InvalidConvexHull;
1155 }
1156
1157 return status() == Status::SimplexCompleted;
1158 }
Index newFacet()
Definition: QuickHull.h:1208
IndexRange pointsOnRidge(const Ridge &R) const
Definition: QuickHull.h:1327
bool areFacetsParallel(const Index if1, const Index if2) const
Definition: QuickHull.h:1273
void deleteFacet(Index f)
Definition: QuickHull.h:1218
Index mergeParallelFacets(const Index if1, const Index if2)
Definition: QuickHull.h:1250
MyPointD Point
Definition: testClone2.cpp:383

References DGtal::QuickHull< TKernel >::above(), DGtal::QuickHull< TKernel >::aboveOrOn(), DGtal::QuickHull< TKernel >::areFacetsNeighbor(), DGtal::QuickHull< TKernel >::areFacetsParallel(), DGtal::QuickHull< TKernel >::assignment, DGtal::QuickHull< TKernel >::checkFacet(), DGtal::QuickHull< TKernel >::checkFacets(), DGtal::QuickHull< TKernel >::checkHull(), DGtal::QuickHull< TKernel >::debug_level, DGtal::QuickHull< TKernel >::deleted_facets, DGtal::QuickHull< TKernel >::deleteFacet(), DGtal::QuickHull< TKernel >::dimension, DGtal::QuickHull< TKernel >::Facet::display(), DGtal::Trace::error(), DGtal::QuickHull< TKernel >::facets, DGtal::H, DGtal::Trace::info(), DGtal::QuickHull< TKernel >::InvalidConvexHull, DGtal::QuickHull< TKernel >::makeFacet(), DGtal::QuickHull< TKernel >::makeNeighbors(), DGtal::QuickHull< TKernel >::mergeParallelFacets(), DGtal::QuickHull< TKernel >::myStatus, DGtal::QuickHull< TKernel >::newFacet(), DGtal::QuickHull< TKernel >::Facet::outside_set, DGtal::QuickHull< TKernel >::pointsOnRidge(), DGtal::QuickHull< TKernel >::processed_points, DGtal::QuickHull< TKernel >::SimplexCompleted, DGtal::QuickHull< TKernel >::status(), DGtal::trace, and DGtal::QuickHull< TKernel >::UNASSIGNED.

Referenced by DGtal::QuickHull< TKernel >::computeFacets().

◆ renumberInfiniteFacets()

template<typename TKernel >
void DGtal::QuickHull< TKernel >::renumberInfiniteFacets ( )
inlineprotected

Determine infinite facets and renumber them so that finite facets come first and infinite facets come after.

Definition at line 926 of file QuickHull.h.

927 {
928 nb_finite_facets = facets.size();
930 if ( ! kernel.hasInfiniteFacets() ) return;
931 IndexRange renumbering( facets.size() );
932 Index i = 0;
933 Index k = facets.size();
934 Index j = 0;
935 for ( auto& l : renumbering ) {
936 if ( ! kernel.isHalfSpaceFacetInfinite( facets[ j ].H ) ) l = i++;
937 else l = --k;
938 j++;
939 }
940 if ( i != k )
941 trace.error() << "[Quickhull::renumberInfiniteFacets]"
942 << " Error renumbering infinite facets "
943 << " up finite=" << i << " low infinite=" << k << std::endl;
944 std::vector< Facet > new_facets( facets.size() );
945 for ( Index f = 0; f < facets.size(); f++ )
946 new_facets[ renumbering[ f ] ].swap( facets[ f ] );
947 facets.swap( new_facets );
948 for ( auto& F : facets ) {
949 for ( auto& N : F.neighbors ) {
950 N = renumbering[ N ];
951 }
952 }
953 // Assign correct number of facets.
956 }
void swap(UnorderedSetByBlock< Key, TSplitter, Hash, KeyEqual, UnorderedMapAllocator > &s1, UnorderedSetByBlock< Key, TSplitter, Hash, KeyEqual, UnorderedMapAllocator > &s2) noexcept

References DGtal::Trace::error(), DGtal::QuickHull< TKernel >::facets, DGtal::QuickHull< TKernel >::kernel, DGtal::QuickHull< TKernel >::nb_finite_facets, DGtal::QuickHull< TKernel >::nb_infinite_facets, DGtal::swap(), and DGtal::trace.

Referenced by DGtal::QuickHull< TKernel >::computeVertices().

◆ selfDisplay()

template<typename TKernel >
void DGtal::QuickHull< TKernel >::selfDisplay ( std::ostream &  out) const
inline

Writes/Displays the object on an output stream.

Parameters
outthe output stream where the object is written.

Definition at line 1519 of file QuickHull.h.

1520 {
1521 out << "[QuickHull " << dimension << "D"
1522 // << " status=" << status()
1523 << " #P=" << nbPoints();
1525 out << " #F=" << nbFacets();
1527 out << " #V=" << nbVertices();
1528 out << "]";
1530 && nbFacets() == 24 ) {
1531 for ( auto f : facets ) f.display( out );
1532 for ( auto v : v2p ) out << points[ v2p[ v ] ] << std::endl;
1533 }
1534 }

References DGtal::QuickHull< TKernel >::AllCompleted, DGtal::QuickHull< TKernel >::dimension, DGtal::QuickHull< TKernel >::facets, DGtal::QuickHull< TKernel >::FacetsCompleted, DGtal::QuickHull< TKernel >::nbFacets(), DGtal::QuickHull< TKernel >::nbPoints(), DGtal::QuickHull< TKernel >::nbVertices(), DGtal::QuickHull< TKernel >::status(), DGtal::QuickHull< TKernel >::v2p, and DGtal::QuickHull< TKernel >::VerticesCompleted.

◆ setInitialSimplex()

template<typename TKernel >
bool DGtal::QuickHull< TKernel >::setInitialSimplex ( const IndexRange full_splx)
inline

Sets the initial full dimensional simplex

Precondition
status() must be Status::InputInitialized
Parameters
full_splxa dimension+1-simplex specified as indices in the vector of input point.
Returns
'true' iff this initial simplex is full dimensional, the object has then the status SimplexCompleted, otherwise returns 'false' and the status is NotFullDimensional.

Definition at line 410 of file QuickHull.h.

411 {
412 if ( status() != Status::InputInitialized ) return false;
413 if ( full_splx.size() != dimension + 1 )
414 {
415 trace.error() << "[QuickHull::setInitialSimplex]"
416 << " not a full dimensional simplex" << std::endl;
418 return false;
419 }
421 for ( Index j = 0; j < dimension; ++j )
422 splx[ j ] = full_splx[ j ];
423 const auto H = kernel.compute( points, splx, full_splx.back() );
424 const auto volume = kernel.volume( H, points[ full_splx.back() ] );
425 if ( volume > 0 )
426 return computeSimplexConfiguration( full_splx );
428 return false;
429 }

References DGtal::QuickHull< TKernel >::computeSimplexConfiguration(), DGtal::QuickHull< TKernel >::dimension, DGtal::Trace::error(), DGtal::H, DGtal::QuickHull< TKernel >::InputInitialized, DGtal::QuickHull< TKernel >::kernel, DGtal::QuickHull< TKernel >::myStatus, DGtal::QuickHull< TKernel >::NotFullDimensional, DGtal::QuickHull< TKernel >::status(), and DGtal::trace.

◆ setInput()

template<typename TKernel >
template<typename InputPoint >
bool DGtal::QuickHull< TKernel >::setInput ( const std::vector< InputPoint > &  input_points,
bool  remove_duplicates = true 
)
inline

Sets the input data for the QuickHull convex hull algorithm, which is a range of points.

Template Parameters
InputPointa model of point that is convertible to Point datatype.
Parameters
[in]input_pointsthe range of input points.
[in]remove_duplicatesshould be set to 'true' if the input data has duplicates.
Returns
'true' if the object is successfully initialized, status must be Status::InputInitialized, 'false' otherwise.
Examples
geometry/tools/exampleLatticeBallQuickHull3D.cpp.

Definition at line 382 of file QuickHull.h.

384 {
385 Clock tic;
386 tic.startClock();
387 clear();
388 timings.clear();
389 kernel.makeInput( points, input2comp, comp2input,
390 input_points, remove_duplicates );
391 timings.push_back( tic.stopClock() );
392 if ( points.size() <= dimension ) {
394 return false;
395 }
397 return true;
398 }
void clear()
Clears the object.
Definition: QuickHull.h:270

References DGtal::QuickHull< TKernel >::clear(), DGtal::QuickHull< TKernel >::comp2input, DGtal::QuickHull< TKernel >::dimension, DGtal::QuickHull< TKernel >::input2comp, DGtal::QuickHull< TKernel >::InputInitialized, DGtal::QuickHull< TKernel >::kernel, DGtal::QuickHull< TKernel >::myStatus, DGtal::QuickHull< TKernel >::NotFullDimensional, tic(), and DGtal::QuickHull< TKernel >::timings.

Referenced by main().

◆ status()

template<typename TKernel >
Status DGtal::QuickHull< TKernel >::status ( ) const
inline

◆ unmakeNeighbors()

template<typename TKernel >
void DGtal::QuickHull< TKernel >::unmakeNeighbors ( const Index  if1,
const Index  if2 
)
inlineprotected

Makes two distinct facets if1 and if2 no more neighbors

Parameters
if1a valid facet index
if2a valid facet index

Definition at line 1238 of file QuickHull.h.

1239 {
1240 facets[ if1 ].subNeighbor( if2 );
1241 facets[ if2 ].subNeighbor( if1 );
1242 }

References DGtal::QuickHull< TKernel >::facets.

Field Documentation

◆ assignment

template<typename TKernel >
std::vector< Index > DGtal::QuickHull< TKernel >::assignment

◆ comp2input

template<typename TKernel >
IndexRange DGtal::QuickHull< TKernel >::comp2input

the injective mapping between the convex hull point range and the input range.

Definition at line 808 of file QuickHull.h.

Referenced by DGtal::QuickHull< TKernel >::clear(), DGtal::QuickHull< TKernel >::memory(), and DGtal::QuickHull< TKernel >::setInput().

◆ debug_level

template<typename TKernel >
int DGtal::QuickHull< TKernel >::debug_level

◆ deleted_facets

template<typename TKernel >
std::set< Index > DGtal::QuickHull< TKernel >::deleted_facets

◆ dimension

template<typename TKernel >
const Size DGtal::QuickHull< TKernel >::dimension = Point::dimension
static

◆ facet_counter

template<typename TKernel >
std::vector< Size > DGtal::QuickHull< TKernel >::facet_counter

Counts the number of facets with a given number of vertices.

Definition at line 828 of file QuickHull.h.

Referenced by DGtal::QuickHull< TKernel >::computeVertices().

◆ facets

template<typename TKernel >
std::vector< Facet > DGtal::QuickHull< TKernel >::facets

◆ input2comp

template<typename TKernel >
IndexRange DGtal::QuickHull< TKernel >::input2comp

the surjective mapping between the input range and the output range used for convex hull computation.

Definition at line 805 of file QuickHull.h.

Referenced by DGtal::QuickHull< TKernel >::clear(), DGtal::QuickHull< TKernel >::memory(), and DGtal::QuickHull< TKernel >::setInput().

◆ kernel

template<typename TKernel >
Kernel DGtal::QuickHull< TKernel >::kernel
mutable

◆ myStatus

template<typename TKernel >
Status DGtal::QuickHull< TKernel >::myStatus
protected

◆ nb_finite_facets

template<typename TKernel >
Size DGtal::QuickHull< TKernel >::nb_finite_facets

Number of finite facets.

Definition at line 822 of file QuickHull.h.

Referenced by DGtal::QuickHull< TKernel >::nbFiniteFacets(), and DGtal::QuickHull< TKernel >::renumberInfiniteFacets().

◆ nb_infinite_facets

template<typename TKernel >
Size DGtal::QuickHull< TKernel >::nb_infinite_facets

Number of infinite facets (!= 0 only for specific kernels)

Definition at line 824 of file QuickHull.h.

Referenced by DGtal::QuickHull< TKernel >::nbInfiniteFacets(), and DGtal::QuickHull< TKernel >::renumberInfiniteFacets().

◆ p2v

template<typename TKernel >
IndexRange DGtal::QuickHull< TKernel >::p2v

◆ points

template<typename TKernel >
std::vector< Point > DGtal::QuickHull< TKernel >::points

the set of points, indexed as in the array.

Definition at line 802 of file QuickHull.h.

◆ processed_points

template<typename TKernel >
IndexRange DGtal::QuickHull< TKernel >::processed_points

◆ timings

template<typename TKernel >
std::vector< double > DGtal::QuickHull< TKernel >::timings

◆ v2p

template<typename TKernel >
IndexRange DGtal::QuickHull< TKernel >::v2p

The documentation for this struct was generated from the following file: