DGtal  1.3.beta
testConvexityHelper.cpp
Go to the documentation of this file.
1 
30 #include <iostream>
32 #include <vector>
33 #include <algorithm>
34 #include "DGtal/base/Common.h"
35 #include "DGtal/kernel/SpaceND.h"
36 #include "DGtal/geometry/volumes/ConvexityHelper.h"
37 #include "DGtal/shapes/SurfaceMesh.h"
38 #include "DGtalCatch.h"
40 
41 using namespace std;
42 using namespace DGtal;
43 
44 
46 // Functions for testing class ConvexityHelper in 2D.
48 
49 SCENARIO( "ConvexityHelper< 2 > unit tests",
50  "[convexity_helper][lattice_polytope][2d]" )
51 {
52  typedef ConvexityHelper< 2 > Helper;
53  typedef Helper::Point Point;
54  typedef ConvexCellComplex< Point > CvxCellComplex;
55  GIVEN( "Given a star { (0,0), (-4,-1), (-3,5), (7,3), (5, -2) } " ) {
56  std::vector<Point> V
57  = { Point(0,0), Point(-4,-1), Point(-3,5), Point(7,3), Point(5, -2) };
58  WHEN( "Computing its lattice polytope" ){
59  const auto P = Helper::computeLatticePolytope( V, false, true );
60  CAPTURE( P );
61  THEN( "The polytope is valid and has 4 non trivial facets" ) {
62  REQUIRE( P.nbHalfSpaces() - 4 == 4 );
63  }
64  THEN( "The polytope is Minkowski summable" ) {
65  REQUIRE( P.canBeSummed() );
66  }
67  THEN( "The polytope contains the input points" ) {
68  REQUIRE( P.isInside( V[ 0 ] ) );
69  REQUIRE( P.isInside( V[ 1 ] ) );
70  REQUIRE( P.isInside( V[ 2 ] ) );
71  REQUIRE( P.isInside( V[ 3 ] ) );
72  REQUIRE( P.isInside( V[ 4 ] ) );
73  }
74  THEN( "The polytope contains 58 points" ) {
75  REQUIRE( P.count() == 58 );
76  }
77  THEN( "The interior of the polytope contains 53 points" ) {
78  REQUIRE( P.countInterior() == 53 );
79  }
80  THEN( "The boundary of the polytope contains 5 points" ) {
81  REQUIRE( P.countBoundary() == 5 );
82  }
83  }
84  }
85  GIVEN( "Given a square with an additional outside vertex " ) {
86  std::vector<Point> V
87  = { Point(-10,-10), Point(10,-10), Point(-10,10), Point(10,10),
88  Point(0,18) };
89  WHEN( "Computing its Delaunay cell complex" ){
90  CvxCellComplex complex;
91  bool ok = Helper::computeDelaunayCellComplex( complex, V, false );
92  CAPTURE( complex );
93  THEN( "The complex has 2 cells, 6 faces, 5 vertices" ) {
94  REQUIRE( ok );
95  REQUIRE( complex.nbCells() == 2 );
96  REQUIRE( complex.nbFaces() == 6 );
97  REQUIRE( complex.nbVertices() == 5 );
98  }
99  THEN( "The faces of cells are finite" ) {
100  bool ok_finite = true;
101  for ( auto c = 0; c < complex.nbCells(); ++c ) {
102  const auto faces = complex.cellFaces( c );
103  for ( auto f : faces )
104  ok_finite = ok_finite && ! complex.isInfinite( complex.faceCell( f ) );
105  }
106  REQUIRE( ok_finite );
107  }
108  THEN( "The opposite of faces of cells are infinite except two" ) {
109  int nb_finite = 0;
110  for ( auto c = 0; c < complex.nbCells(); ++c ) {
111  const auto faces = complex.cellFaces( c );
112  for ( auto f : faces ) {
113  const auto opp_f = complex.opposite( f );
114  nb_finite += complex.isInfinite( complex.faceCell( opp_f ) ) ? 0 : 1;
115  }
116  }
117  REQUIRE( nb_finite == 2 );
118  }
119  }
120  }
121  GIVEN( "Given a degenerated polytope { (0,0), (3,-1), (9,-3), (-6,2) } " ) {
122  std::vector<Point> V
123  = { Point(0,0), Point(3,-1), Point(9,-3), Point(-6,2) };
124  WHEN( "Computing its lattice polytope" ){
125  const auto P = Helper::computeLatticePolytope( V, false, true );
126  CAPTURE( P );
127  THEN( "The polytope is valid and has 2 non trivial facets" ) {
128  REQUIRE( P.nbHalfSpaces() - 4 == 2 );
129  }
130  THEN( "The polytope contains 6 points" ) {
131  REQUIRE( P.count() == 6 );
132  }
133  THEN( "The polytope contains no interior points" ) {
134  REQUIRE( P.countInterior() == 0 );
135  }
136  }
137  }
138  GIVEN( "Given a degenerated simplex { (4,0), (7,2), (-5,-6) } " ) {
139  std::vector<Point> V
140  = { Point(4,0), Point(7,2), Point(-5,-6) };
141  WHEN( "Computing its lattice polytope" ){
142  const auto P = Helper::computeLatticePolytope( V, false, true );
143  CAPTURE( P );
144  THEN( "The polytope is valid and has 2 non trivial facets" ) {
145  REQUIRE( P.nbHalfSpaces() - 4 == 2 );
146  }
147  THEN( "The polytope contains 5 points" ) {
148  REQUIRE( P.count() == 5 );
149  }
150  THEN( "The polytope contains no interior points" ) {
151  REQUIRE( P.countInterior() == 0 );
152  }
153  }
154  }
155 }
156 
158 // Functions for testing class ConvexityHelper in 3D.
160 
161 SCENARIO( "ConvexityHelper< 3 > unit tests",
162  "[convexity_helper][3d]" )
163 {
164  typedef ConvexityHelper< 3 > Helper;
165  typedef Helper::Point Point;
169  typedef PolygonalSurface< Point > LatticePolySurf;
170  typedef ConvexCellComplex< Point > CvxCellComplex;
171  GIVEN( "Given an octahedron star { (0,0,0), (-2,0,0), (2,0,0), (0,-2,0), (0,2,0), (0,0,-2), (0,0,2) } " ) {
172  std::vector<Point> V
173  = { Point(0,0,0), Point(-2,0,0), Point(2,0,0), Point(0,-2,0), Point(0,2,0),
174  Point(0,0,-2), Point(0,0,2) };
175  WHEN( "Computing its lattice polytope" ){
176  const auto P = Helper::computeLatticePolytope( V, false, true );
177  CAPTURE( P );
178  THEN( "The polytope is valid and has 8 non trivial facets plus 12 edge constraints" ) {
179  REQUIRE( P.nbHalfSpaces() - 6 == 20 );
180  }
181  THEN( "The polytope is Minkowski summable" ) {
182  REQUIRE( P.canBeSummed() );
183  }
184  THEN( "The polytope contains the input points" ) {
185  REQUIRE( P.isInside( V[ 0 ] ) );
186  REQUIRE( P.isInside( V[ 1 ] ) );
187  REQUIRE( P.isInside( V[ 2 ] ) );
188  REQUIRE( P.isInside( V[ 3 ] ) );
189  REQUIRE( P.isInside( V[ 4 ] ) );
190  REQUIRE( P.isInside( V[ 5 ] ) );
191  REQUIRE( P.isInside( V[ 6 ] ) );
192  }
193  THEN( "The polytope contains 25 points" ) {
194  REQUIRE( P.count() == 25 );
195  }
196  THEN( "The interior of the polytope contains 7 points" ) {
197  REQUIRE( P.countInterior() == 7 );
198  }
199  THEN( "The boundary of the polytope contains 18 points" ) {
200  REQUIRE( P.countBoundary() == 18 );
201  }
202  }
203  WHEN( "Computing the boundary of its convex hull as a SurfaceMesh" ){
204  SMesh smesh;
205  bool ok = Helper::computeConvexHullBoundary( smesh, V, false );
206  CAPTURE( smesh );
207  THEN( "The surface mesh is valid and has 6 vertices, 12 edges and 8 faces" ) {
208  REQUIRE( ok );
209  REQUIRE( smesh.nbVertices() == 6 );
210  REQUIRE( smesh.nbEdges() == 12 );
211  REQUIRE( smesh.nbFaces() == 8 );
212  }
213  THEN( "The surface mesh has the topology of a sphere" ) {
214  REQUIRE( smesh.Euler() == 2 );
215  REQUIRE( smesh.computeManifoldBoundaryEdges().size() == 0 );
216  REQUIRE( smesh.computeNonManifoldEdges().size() == 0 );
217  }
218  }
219  WHEN( "Computing the boundary of its convex hull as a lattice PolygonalSurface" ){
220  LatticePolySurf lpsurf;
221  bool ok = Helper::computeConvexHullBoundary( lpsurf, V, false );
222  CAPTURE( lpsurf );
223  THEN( "The polygonal surface is valid and has 6 vertices, 12 edges and 8 faces" ) {
224  REQUIRE( ok );
225  REQUIRE( lpsurf.nbVertices() == 6 );
226  REQUIRE( lpsurf.nbEdges() == 12 );
227  REQUIRE( lpsurf.nbFaces() == 8 );
228  REQUIRE( lpsurf.nbArcs() == 24 );
229  }
230  THEN( "The polygonal surface has the topology of a sphere and no boundary" ) {
231  REQUIRE( lpsurf.Euler() == 2 );
232  REQUIRE( lpsurf.allBoundaryArcs().size() == 0 );
233  REQUIRE( lpsurf.allBoundaryVertices().size() == 0 );
234  }
235  }
236  WHEN( "Computing its convex hull as a ConvexCellComplex" ){
237  CvxCellComplex complex;
238  bool ok = Helper::computeConvexHullCellComplex( complex, V, false );
239  CAPTURE( complex );
240  THEN( "The convex cell complex is valid and has 6 vertices, 8 faces and 1 finite cell" ) {
241  REQUIRE( ok );
242  REQUIRE( complex.nbVertices() == 6 );
243  REQUIRE( complex.nbFaces() == 8 );
244  REQUIRE( complex.nbCells() == 1 );
245  }
246  }
247  }
248  GIVEN( "Given a cube with an additional outside vertex " ) {
249  std::vector<Point> V
250  = { Point(-10,-10,-10), Point(10,-10,-10), Point(-10,10,-10), Point(10,10,-10),
251  Point(-10,-10,10), Point(10,-10,10), Point(-10,10,10), Point(10,10,10),
252  Point(0,0,18) };
253  WHEN( "Computing its Delaunay cell complex" ){
254  CvxCellComplex complex;
255  bool ok = Helper::computeDelaunayCellComplex( complex, V, false );
256  CAPTURE( complex );
257  THEN( "The complex has 2 cells, 10 faces, 9 vertices" ) {
258  REQUIRE( ok );
259  REQUIRE( complex.nbCells() == 2 );
260  REQUIRE( complex.nbFaces() == 10 );
261  REQUIRE( complex.nbVertices() == 9 );
262  }
263  THEN( "The faces of cells are finite" ) {
264  bool ok_finite = true;
265  for ( auto c = 0; c < complex.nbCells(); ++c ) {
266  const auto faces = complex.cellFaces( c );
267  for ( auto f : faces )
268  ok_finite = ok_finite && ! complex.isInfinite( complex.faceCell( f ) );
269  }
270  REQUIRE( ok_finite );
271  }
272  THEN( "The opposite of faces of cells are infinite except two" ) {
273  int nb_finite = 0;
274  for ( auto c = 0; c < complex.nbCells(); ++c ) {
275  const auto faces = complex.cellFaces( c );
276  for ( auto f : faces ) {
277  const auto opp_f = complex.opposite( f );
278  nb_finite += complex.isInfinite( complex.faceCell( opp_f ) ) ? 0 : 1;
279  }
280  }
281  REQUIRE( nb_finite == 2 );
282  }
283  }
284  }
285  GIVEN( "Given a degenerated 1d polytope { (0,0,1), (3,-1,2), (9,-3,4), (-6,2,-1) } " ) {
286  std::vector<Point> V
287  = { Point(0,0,1), Point(3,-1,2), Point(9,-3,4), Point(-6,2,-1) };
288  WHEN( "Computing its lattice polytope" ){
289  const auto P = Helper::computeLatticePolytope( V, false, true );
290  CAPTURE( P );
291  THEN( "The polytope is valid and has 6 non trivial facets" ) {
292  REQUIRE( P.nbHalfSpaces() - 6 == 6 );
293  }
294  THEN( "The polytope contains 6 points" ) {
295  REQUIRE( P.count() == 6 );
296  }
297  THEN( "The polytope contains no interior points" ) {
298  REQUIRE( P.countInterior() == 0 );
299  }
300  }
301  }
302  GIVEN( "Given a degenerated 1d simplex { (1,0,-1), Point(4,-1,-2), Point(10,-3,-4) } " ) {
303  std::vector<Point> V
304  = { Point(1,0,-1), Point(4,-1,-2), Point(10,-3,-4) };
305  WHEN( "Computing its lattice polytope" ){
306  const auto P = Helper::computeLatticePolytope( V, false, true );
307  CAPTURE( P );
308  THEN( "The polytope is valid and has 6 non trivial facets" ) {
309  REQUIRE( P.nbHalfSpaces() - 6 == 6 );
310  }
311  THEN( "The polytope contains 4 points" ) {
312  REQUIRE( P.count() == 4 );
313  }
314  THEN( "The polytope contains no interior points" ) {
315  REQUIRE( P.countInterior() == 0 );
316  }
317  }
318  }
319  GIVEN( "Given a degenerated 2d polytope { (2,1,0), (1,0,1), (1,2,1), (0,1,2), (0,3,2) } " ) {
320  std::vector<Point> V
321  = { Point(2,1,0), Point(1,0,1), Point(1,2,1), Point(0,1,2), Point(0,3,2) };
322  WHEN( "Computing its lattice polytope" ){
323  const auto P = Helper::computeLatticePolytope( V, false, true );
324  CAPTURE( P );
325  THEN( "The polytope is valid and has more than 6 non trivial facets" ) {
326  REQUIRE( P.nbHalfSpaces() - 6 == 6 );
327  }
328  THEN( "The polytope contains 7 points" ) {
329  REQUIRE( P.count() == 7 );
330  }
331  THEN( "The polytope contains no interior points" ) {
332  REQUIRE( P.countInterior() == 0 );
333  }
334  }
335  }
336  GIVEN( "Given a degenerated 2d simplex { (2,1,0), (1,0,1), (1,5,1), (0,3,2) } " ) {
337  std::vector<Point> V
338  = { Point(2,1,0), Point(1,0,1), Point(1,5,1), Point(0,3,2) };
339  WHEN( "Computing its lattice polytope" ){
340  const auto P = Helper::computeLatticePolytope( V, false, true );
341  CAPTURE( P );
342  THEN( "The polytope is valid and has more than 6 non trivial facets" ) {
343  REQUIRE( P.nbHalfSpaces() - 6 == 6 );
344  }
345  THEN( "The polytope contains 8 points" ) {
346  REQUIRE( P.count() == 8 );
347  }
348  THEN( "The polytope contains no interior points" ) {
349  REQUIRE( P.countInterior() == 0 );
350  }
351  }
352  }
353 }
DGtal::SurfaceMesh::Euler
long Euler() const
Definition: SurfaceMesh.h:294
DGtal::SurfaceMesh
Aim: Represents an embedded mesh as faces and a list of vertices. Vertices may be shared among faces ...
Definition: SurfaceMesh.h:91
DGtal::SurfaceMesh::computeManifoldBoundaryEdges
Edges computeManifoldBoundaryEdges() const
REQUIRE
REQUIRE(domain.isInside(aPoint))
DGtal::SurfaceMesh::nbVertices
Size nbVertices() const
Definition: SurfaceMesh.h:280
DGtal::ConvexCellComplex
Aim: represents a d-dimensional complex in a d-dimensional space with the following properties and re...
Definition: ConvexCellComplex.h:85
CAPTURE
CAPTURE(thicknessHV)
DGtal::SurfaceMesh::computeNonManifoldEdges
Edges computeNonManifoldEdges() const
DGtal::SurfaceMesh::nbEdges
Size nbEdges() const
Definition: SurfaceMesh.h:284
DGtal
DGtal is the top-level namespace which contains all DGtal functions and types.
RealVector
Space::RealVector RealVector
Definition: fullConvexitySphereGeodesics.cpp:113
DGtal::ConvexityHelper
Aim: Provides a set of functions to facilitate the computation of convex hulls and polytopes,...
Definition: ConvexityHelper.h:164
SCENARIO
SCENARIO("ConvexityHelper< 2 > unit tests", "[convexity_helper][lattice_polytope][2d]")
Definition: testConvexityHelper.cpp:49
SMesh
SurfaceMesh< RealPoint, RealVector > SMesh
Definition: fullConvexitySphereGeodesics.cpp:115
GIVEN
GIVEN("A cubical complex with random 3-cells")
Definition: testCubicalComplex.cpp:70
Point
MyPointD Point
Definition: testClone2.cpp:383
DGtal::PolygonalSurface
Aim: Represents a polygon mesh, i.e. a 2-dimensional combinatorial surface whose faces are (topologic...
Definition: PolygonalSurface.h:86
RealPoint
Z2i::RealPoint RealPoint
Definition: testAstroid2D.cpp:46
DGtal::SurfaceMesh::nbFaces
Size nbFaces() const
Definition: SurfaceMesh.h:288