CGAL::convex_hull_incremental_3

Definition

The function convex_hull_incremental_3 computes the convex hull polyhedron from a set of given three-dimensional points.

#include <CGAL/convex_hull_incremental_3.h>

template <class InputIterator, class Polyhedron>
void
convex_hull_incremental_3 ( InputIterator first,
InputIterator beyond,
Polyhedron& P,
bool test_correctness = false)
computes the convex hull polyhedron of the points in the range [first,beyond) and assigns it to P. If test_correctness is set to true, the tests described in [MNS+96] are used to determine the correctness of the resulting polyhedron.

Requirements

  1. Polyhedron must provide a type Polyhedron::Traits that defines the following types
  2. InputIterator::value_type must be the same as Polyhedron::Traits::Point

See Also

CGAL::convex_hull_3
CGAL::convex_hull_2

Implementation

This function uses the d-dimensional convex hull incremental construction algorithm [CMS93] with d fixed to 3. The algorithm requires O(n2) time in the worst case and O(n log n) expected time.

See Also

CGAL::Convex_hull_d<R>

Example

The following example computes the convex hull of a set of 250 random points chosen uniformly in a sphere of radius 100.

File: examples/Convex_hull_3/incremental_hull_3.cpp
#include <CGAL/Homogeneous.h>
#include <CGAL/point_generators_3.h>
#include <CGAL/algorithm.h>
#include <CGAL/Polyhedron_3.h>
#include <CGAL/convex_hull_incremental_3.h>
#include <vector>

#ifdef CGAL_USE_GMP
#include <CGAL/Gmpz.h>
typedef CGAL::Gmpz RT;
#else
#include <CGAL/MP_Float.h>
typedef CGAL::MP_Float RT;
#endif


typedef CGAL::Homogeneous<RT>                  K;
typedef K::Point_3                             Point_3;
typedef CGAL::Polyhedron_3< K>                 Polyhedron;
typedef CGAL::Creator_uniform_3<int, Point_3>  Creator;

int main()
{
  CGAL::Random_points_in_sphere_3<Point_3, Creator> gen(100.0);

  std::vector<Point_3> V;
  // generate 250 points randomly on a sphere of radius 100.0 and copy
  // them to a vector
  CGAL::copy_n( gen, 250, std::back_inserter(V) );

  Polyhedron P; // define polyhedron to hold convex hull

  // compute convex hull
  CGAL::convex_hull_incremental_3( V.begin(), V.end(), P, true);


  return 0;
}