DGtal  1.5.beta
exampleRationalBallQuickHull3D.cpp File Reference
#include "DGtal/base/Common.h"
#include "DGtal/geometry/tools/QuickHull.h"
#include "DGtal/shapes/SurfaceMesh.h"
#include "DGtal/io/writers/SurfaceMeshWriter.h"
Include dependency graph for exampleRationalBallQuickHull3D.cpp:

Go to the source code of this file.

Functions

double rand01 ()
 [QuickHull3D-Includes] More...
 
int main (int argc, char *argv[])
 

Detailed Description

This program is free software: you can redistribute it and/or modify it under the terms of the GNU Lesser General Public License as published by the Free Software Foundation, either version 3 of the License, or (at your option) any later version.

This program is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License for more details.

You should have received a copy of the GNU General Public License along with this program. If not, see http://www.gnu.org/licenses/.

Author
Jacques-Olivier Lachaud (jacqu.nosp@m.es-o.nosp@m.livie.nosp@m.r.la.nosp@m.chaud.nosp@m.@uni.nosp@m.v-sav.nosp@m.oie..nosp@m.fr ) Laboratory of Mathematics (CNRS, UMR 5127), University of Savoie, France
Date
2021/01/01

An example file named exampleRationalBallQuickHull3D.

This file is part of the DGtal library.

Definition in file exampleRationalBallQuickHull3D.cpp.

Function Documentation

◆ main()

int main ( int  argc,
char *  argv[] 
)

[QuickHull3D-Typedefs]

[QuickHull3D-Typedefs]

[QuickHull3D-Computation]

[QuickHull3D-Computation] [QuickHull3D-Timings]

[QuickHull3D-Timings]

[QuickHull3D-BuildMesh]

[QuickHull3D-BuildMesh]

[QuickHull3D-OutputMesh]

[QuickHull3D-OutputMesh]

Definition at line 63 of file exampleRationalBallQuickHull3D.cpp.

64 {
65  int nb = argc > 1 ? atoi( argv[ 1 ] ) : 100; // nb points
66  double R = argc > 2 ? atof( argv[ 2 ] ) : 10.0; // radius of ball
67  double precision = argc > 3 ? atof( argv[ 3 ] ) : 1024.0; // precision
68  // (0) typedefs
73  // (1) create range of random points in ball
74  std::vector< RealPoint > V;
75  const double R2 = R * R;
76  for ( int i = 0; i < nb; ) {
77  RealPoint p( rand01()*2.0*R - R, rand01()*2.0*R - R, rand01()*2.0*R - R );
78  if ( p.squaredNorm() < R2 ) { V.push_back( p ); i++; }
79  }
80  // (2) compute convex hull
82  Kernel3D kernel( precision );
83  QuickHull3D hull( kernel );
84  hull.setInput( V );
85  hull.computeConvexHull();
86  std::cout << "#points=" << hull.nbPoints()
87  << " #vertices=" << hull.nbVertices()
88  << " #facets=" << hull.nbFacets() << std::endl;
91  double total_time = 0;
92  std::for_each( hull.timings.cbegin(), hull.timings.cend(),
93  [&total_time] ( double t ) { total_time += t; } );
94  std::cout << "purge duplicates= " << round(hull.timings[ 0 ]) << " ms." << std::endl;
95  std::cout << "init simplex = " << round(hull.timings[ 1 ]) << " ms." << std::endl;
96  std::cout << "quickhull core = " << round(hull.timings[ 2 ]) << " ms." << std::endl;
97  std::cout << "compute vertices= " << round(hull.timings[ 3 ]) << " ms." << std::endl;
98  std::cout << "total time = " << round(total_time) << " ms." << std::endl;
100  // (3) build mesh
102  std::vector< RealPoint > positions;
103  hull.getVertexPositions( positions );
104  std::vector< std::vector< std::size_t > > facets;
105  hull.getFacetVertices( facets );
107  SMesh mesh( positions.cbegin(), positions.cend(), facets.cbegin(), facets.cend() );
109  // (4) output result as OBJ file
111  std::ofstream out( "qhull.obj" );
113  out.close();
115  return 0;
116 }
Aim: Implements basic operations that will be used in Point and Vector classes.
Definition: PointVector.h:593
double rand01()
[QuickHull3D-Includes]
SurfaceMesh< RealPoint, RealVector > SMesh
Aim: a geometric kernel to compute the convex hull of digital points with integer-only arithmetic.
Aim: a geometric kernel to compute the convex hull of floating points with integer-only arithmetic....
Aim: Implements the quickhull algorithm by Barber et al. , a famous arbitrary dimensional convex hull...
Definition: QuickHull.h:140
static bool writeOBJ(std::ostream &output, const SurfaceMesh &smesh)
Aim: Represents an embedded mesh as faces and a list of vertices. Vertices may be shared among faces ...
Definition: SurfaceMesh.h:92

References DGtal::QuickHull< TKernel >::computeConvexHull(), DGtal::QuickHull< TKernel >::getFacetVertices(), DGtal::QuickHull< TKernel >::getVertexPositions(), DGtal::QuickHull< TKernel >::nbFacets(), DGtal::QuickHull< TKernel >::nbPoints(), DGtal::QuickHull< TKernel >::nbVertices(), rand01(), DGtal::QuickHull< TKernel >::setInput(), DGtal::QuickHull< TKernel >::timings, and DGtal::SurfaceMeshWriter< TRealPoint, TRealVector >::writeOBJ().

◆ rand01()

double rand01 ( )

[QuickHull3D-Includes]

[QuickHull3D-Includes]

Definition at line 60 of file exampleRationalBallQuickHull3D.cpp.

60 { return (double) rand() / (double) RAND_MAX; }

Referenced by main().