songsenand
/
cgal
mirror of https://github.com/CGAL/cgal
1
0
Fork 0
Code Issues Packages Projects Releases Wiki Activity

use boost::optional for the computation (or not) of an interior point

This commit is contained in:
Jocelyn MEYRON 2014-11-11 15:51:20 +01:00
parent 1bb886e9dc
commit b452d5a7e1
6 changed files with 120 additions and 100 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

4
Convex_hull_3/doc/Convex_hull_3/CGAL/Convex_hull_3/dual/halfspace_intersection_3.h
View File

@ -4,7 +4,7 @@ namespace CGAL {
\ingroup PkgConvexHull3Functions
\brief computes robustly the intersection of the halfspaces defined by the planes contained in the range [`begin`, `end`) without constructing the dual points. The result is stored in the polyhedron `P`.
`origin` is a point strictly inside the polyhedron.
If `origin` is given then it must be a point strictly inside the polyhedron. If an interior point is not given then it is computed using a linear program and thus is slower than the first approach.
This version does not compute the dual points, but instead it uses a traits class for handling predicates for dual points without computing them.
\attention Halfspaces are considered as lower halfspaces that is to say if the plane's equation is \f$ a\, x +b\, y +c\, z + d = 0 \f$ then the corresponding halfspace is defined by \f$ a\, x +b\, y +c\, z + d \le 0 \f$ .
@ -22,6 +22,6 @@ This version does not compute the dual points, but instead it uses a traits clas
template <class PlaneIterator, class Polyhedron>
void halfspace_intersection_3 (PlaneIterator begin, PlaneIterator end,
Polyhedron &P,
typename Polyhedron::Vertex::Point_3 const& origin);
Point_3 origin);
} /* namespace CGAL */

6
Convex_hull_3/doc/Convex_hull_3/CGAL/Convex_hull_3/dual/halfspace_intersection_with_constructions_3.h
View File

@ -1,10 +1,10 @@
namespace CGAL {
/*!
\ingroup PkgConvexHull3Functions
\Ingroup PkgConvexHull3Functions
\brief computes the intersection of the halfspaces defined by the planes contained in the range [`begin`, `end`). The result is stored in the polyhedron `P`.
`origin` is a point strictly inside the polyhedron.
If `origin` is given then it must be a point strictly inside the polyhedron. If an interior point is not given then it is computed using a linear program and thus is slower than the first approach.
This version constructs explicitly the dual points using the convex hull algorithm parametrized with the given traits class.
\attention Halfspaces are considered as lower halfspaces that is to say if the plane's equation is \f$ a\, x +b\, y +c\, z + d = 0 \f$ then the corresponding halfspace is defined by \f$ a\, x +b\, y +c\, z + d \le 0 \f$ .
@ -24,7 +24,7 @@ template <class PlaneIterator, class Polyhedron, class Traits>
void halfspace_intersection_with_constructions_3(PlaneIterator pbegin,
PlaneIterator pend,
Polyhedron &P,
typename Polyhedron::Vertex::Point_3 const& origin,
Point_3 origin,
const Traits & ch_traits = Default_traits);
} /* namespace CGAL */

6
Convex_hull_3/doc/Convex_hull_3/Convex_hull_3.txt
View File

@ -67,11 +67,15 @@ point set or not. This function is used in postcondition testing for
\subsection Convex_hull_3HalfspaceIntersection Halfspace Intersection
The functions `halfspace_intersection_3()` and
`halfspace_intersection_with_constructions_3()`
uses the convex hull computation and the duality to compute the intersection
uses the convex hull algorithm and the duality to compute the intersection
of a list of halfspaces. The first version does not explicitly computes
the dual points: the traits class handles this issue. The second one constructs
these points and hence is less robust but the computation is faster.
In order to compute, the intersection an interior point is needed. It could be
either given by the user or computed using linear programming. Notice that the
second approach is slower due to the resolution of a linear program.
The following program computes the intersection of halfspaces defined by
tangent planes to a sphere:

2
Convex_hull_3/examples/Convex_hull_3/halfspace_intersection_3.cpp
View File

@ -38,7 +38,7 @@ int main (void) {
planes.end(),
P,
Point(0, 0, 0));
return 0;
}

91
Convex_hull_3/include/CGAL/Convex_hull_3/dual/halfspace_intersection_3.h
View File

@ -257,37 +257,13 @@ namespace CGAL
} // namespace internal
} // namespace Convex_hull_3
// Compute the intersection of halfspaces by computing a point
// inside the polyhedron (using linear programming)
template <class PlaneIterator, class Polyhedron>
void halfspace_intersection_without_origin_3 (PlaneIterator begin, PlaneIterator end,
Polyhedron &P) {
// Types
typedef typename Kernel_traits<typename Polyhedron::Vertex::Point_3>::Kernel K;
typedef typename Polyhedron::Vertex::Point_3 Point_3;
// choose exact integral type
#ifdef CGAL_USE_GMP
typedef CGAL::Gmpq ET;
#else
typedef CGAL::MP_Float ET;
#endif
// find a point inside the intersection
typedef Interior_polyhedron_3<K, ET> Interior_polyhedron;
Interior_polyhedron interior;
bool res = interior.find(begin, end);
CGAL_assertion_msg(res, "halfspace_intersection_without_origin_3: problem when determing a point inside");
Point_3 origin = interior.inside_point();
// compute the intersection
halfspace_intersection_3(begin, end, P, origin);
}
// Compute the intersection of halfspaces with the origin given by the user
// Compute the intersection of halfspaces.
// If the user gives an origin then the function used it, otherwise, it is
// computed using a linear program.
template <class PlaneIterator, class Polyhedron>
void halfspace_intersection_3 (PlaneIterator begin, PlaneIterator end,
Polyhedron &P,
typename Polyhedron::Vertex::Point_3 const& origin) {
boost::optional<typename Polyhedron::Vertex::Point_3> const& origin) {
// Checks whether the intersection if a polyhedron
CGAL_assertion_msg(Convex_hull_3::internal::is_intersection_dim_3(begin, end), "halfspace_intersection_3: intersection not a polyhedron");
@ -296,18 +272,61 @@ namespace CGAL
typedef Convex_hull_3::Convex_hull_traits_dual_3<K> Hull_traits_dual_3;
typedef Polyhedron_3<Hull_traits_dual_3> Polyhedron_dual_3;
typedef Convex_hull_3::internal::Build_primal_polyhedron<K, Polyhedron_dual_3, Polyhedron> Builder;
typedef typename Polyhedron::Vertex::Point_3 Point_3;
Hull_traits_dual_3 dual_traits(origin);
if (origin) {
Point_3 p_origin = boost::get(origin);
Hull_traits_dual_3 dual_traits(p_origin);
Polyhedron_dual_3 dual_convex_hull;
CGAL::convex_hull_3(begin, end, dual_convex_hull, dual_traits);
Builder build_primal(dual_convex_hull, origin);
P.delegate(build_primal);
Polyhedron_dual_3 dual_convex_hull;
CGAL::convex_hull_3(begin, end, dual_convex_hull, dual_traits);
Builder build_primal(dual_convex_hull, p_origin);
P.delegate(build_primal);
// Posterior check if the origin is inside the computed polyhedron
Polyhedron Q(P);
CGAL_assertion_msg(!Convex_hull_3::internal::point_inside_convex_polyhedron(Q, origin), "halfspace_intersection_3: origin not in the polyhedron");
// Posterior check if the origin is inside the computed polyhedron
Polyhedron Q(P);
CGAL_assertion_msg(!Convex_hull_3::internal::point_inside_convex_polyhedron(Q, p_origin), "halfspace_intersection_3: origin not in the polyhedron");
} else {
// choose exact integral type
#ifdef CGAL_USE_GMP
typedef CGAL::Gmpq ET;
#else
typedef CGAL::MP_Float ET;
#endif
// find a point inside the intersection
typedef Interior_polyhedron_3<K, ET> Interior_polyhedron;
Interior_polyhedron interior;
bool res = interior.find(begin, end);
CGAL_assertion_msg(res, "halfspace_intersection_without_origin_3: problem when determing a point inside");
Point_3 origin = interior.inside_point();
Hull_traits_dual_3 dual_traits(origin);
Polyhedron_dual_3 dual_convex_hull;
CGAL::convex_hull_3(begin, end, dual_convex_hull, dual_traits);
Builder build_primal(dual_convex_hull, origin);
P.delegate(build_primal);
// Posterior check if the origin is inside the computed polyhedron
Polyhedron Q(P);
CGAL_assertion_msg(!Convex_hull_3::internal::point_inside_convex_polyhedron(Q, origin), "halfspace_intersection_3: origin not in the polyhedron");
}
}
// Compute the intersection of halfspaces by computing a point inside.
template <class PlaneIterator, class Polyhedron>
void halfspace_intersection_3 (PlaneIterator begin, PlaneIterator end,
Polyhedron &P) {
halfspace_intersection_3(begin, end , P, boost::none);
}
// Compute the intersection of halfspaces (an interior point is given.)
template <class PlaneIterator, class Polyhedron>
void halfspace_intersection_3 (PlaneIterator begin, PlaneIterator end,
Polyhedron &P,
typename Polyhedron::Vertex::Point_3 const& origin) {
halfspace_intersection_3(begin, end , P, boost::optional<typename Polyhedron::Vertex::Point_3>(origin));
}
} // namespace CGAL
#endif // CGAL_HALFSPACE_INTERSECTION_3_H

111
Convex_hull_3/include/CGAL/Convex_hull_3/dual/halfspace_intersection_with_constructions_3.h
View File

@ -116,14 +116,32 @@ namespace CGAL
void halfspace_intersection_with_constructions_3(PlaneIterator pbegin,
PlaneIterator pend,
Polyhedron &P,
typename Polyhedron::Vertex::Point_3 const& origin,
const Traits & ch_traits)
{
boost::optional<typename Polyhedron::Vertex::Point_3> const& origin,
const Traits & ch_traits) {
typedef typename Kernel_traits<typename Polyhedron::Vertex::Point_3>::Kernel K;
typedef typename K::Point_3 Point;
typedef typename K::Plane_3 Plane;
typedef typename CGAL::internal::Build_dual_polyhedron<Polyhedron> Builder;
Point p_origin;
if (origin) {
p_origin = boost::get(origin);
} else {
// choose exact integral type
#ifdef CGAL_USE_GMP
typedef CGAL::Gmpq ET;
#else
typedef CGAL::MP_Float ET;
#endif
// find a point inside the intersection
typedef Interior_polyhedron_3<K, ET> Interior_polyhedron;
Interior_polyhedron interior;
bool res = interior.find(pbegin, pend);
CGAL_assertion_msg(res, "halfspace_intersection_with_constructions_3: problem when determing a point inside");
p_origin = interior.inside_point();
}
// construct dual points to apply the convex hull
std::vector<Point> dual_points;
for (PlaneIterator p = pbegin; p != pend; ++p) {
@ -131,54 +149,39 @@ namespace CGAL
Plane translated_p(p->a(),
p->b(),
p->c(),
p->d() + origin.x() * p->a() + origin.y() * p->b() + origin.z() * p->c());
p->d() + p_origin.x() * p->a() + p_origin.y() * p->b() + p_origin.z() * p->c());
dual_points.push_back(CGAL::ORIGIN + translated_p.orthogonal_vector () / (-translated_p.d()));
}
Polyhedron ch;
CGAL::convex_hull_3(dual_points.begin(), dual_points.end(), ch, ch_traits);
Builder build_dual (ch, origin);
Builder build_dual (ch, p_origin);
P.delegate(build_dual);
}
// Compute the intersection of halfspaces by constructing explicitly
// the dual points with the traits class for convex_hull_3 given
// as an argument.
// The point inside the polyhedron is computed using linear programming.
// An interior point is given.
template <class PlaneIterator, class Polyhedron, class Traits>
void halfspace_intersection_with_constructions_without_origin_3 (PlaneIterator begin, PlaneIterator end,
Polyhedron &P,
const Traits & ch_traits) {
// Types
typedef typename Kernel_traits<typename Polyhedron::Vertex::Point_3>::Kernel K;
typedef typename Polyhedron::Vertex::Point_3 Point_3;
// choose exact integral type
#ifdef CGAL_USE_GMP
typedef CGAL::Gmpq ET;
#else
typedef CGAL::MP_Float ET;
#endif
// find a point inside the intersection
typedef Interior_polyhedron_3<K, ET> Interior_polyhedron;
Interior_polyhedron interior;
bool res = interior.find(begin, end);
CGAL_assertion_msg(res, "halfspace_intersection_with_constructions_without_origin_3: problem when determing a point inside");
Point_3 origin = interior.inside_point();
// compute the intersection
halfspace_intersection_with_constructions_3(begin, end, P, origin, ch_traits);
}
// Compute the intersection of halfspaces by constructing explicitly
// the dual points with the default traits class for convex_hull_3
template <class PlaneIterator, class Polyhedron>
void halfspace_intersection_with_constructions_3(PlaneIterator pbegin,
PlaneIterator pend,
Polyhedron &P,
typename Polyhedron::Vertex::Point_3 const& origin)
{
typename Polyhedron::Vertex::Point_3 const& origin,
const Traits & ch_traits) {
halfspace_intersection_with_constructions_3(pbegin, pend, P,
boost::optional<typename Polyhedron::Vertex::Point_3>(origin),
ch_traits);
}
// Compute the intersection of halfspaces by constructing explicitly
// the dual points with the default traits class for convex_hull_3.
template <class PlaneIterator, class Polyhedron>
void halfspace_intersection_with_constructions_3 (PlaneIterator pbegin,
PlaneIterator pend,
Polyhedron &P,
boost::optional<typename Polyhedron::Vertex::Point_3> const& origin) {
typedef typename Kernel_traits<typename Polyhedron::Vertex::Point_3>::Kernel K;
typedef typename K::Point_3 Point_3;
typedef typename internal::Convex_hull_3::Default_traits_for_Chull_3<Point_3>::type Traits;
@ -188,30 +191,24 @@ namespace CGAL
// Compute the intersection of halfspaces by constructing explicitly
// the dual points with the default traits class for convex_hull_3.
// The point inside the polyhedron is computed using linear programming.
// An interior point is given.
template <class PlaneIterator, class Polyhedron>
void halfspace_intersection_with_constructions_without_origin_3 (PlaneIterator begin, PlaneIterator end,
Polyhedron &P) {
// Types
typedef typename Kernel_traits<typename Polyhedron::Vertex::Point_3>::Kernel K;
typedef typename Polyhedron::Vertex::Point_3 Point_3;
typedef typename internal::Convex_hull_3::Default_traits_for_Chull_3<Point_3>::type Traits;
void halfspace_intersection_with_constructions_3 (PlaneIterator pbegin,
PlaneIterator pend,
Polyhedron &P,
typename Polyhedron::Vertex::Point_3 const& origin) {
halfspace_intersection_with_constructions_3(pbegin, pend, P,
boost::optional<typename Polyhedron::Vertex::Point_3>(origin));
}
// choose exact integral type
#ifdef CGAL_USE_GMP
typedef CGAL::Gmpq ET;
#else
typedef CGAL::MP_Float ET;
#endif
// find a point inside the intersection
typedef Interior_polyhedron_3<K, ET> Interior_polyhedron;
Interior_polyhedron interior;
bool res = interior.find(begin, end);
CGAL_assertion_msg(res, "halfspace_intersection_with_constructions_without_origin_3: problem when determing an point inside");
Point_3 origin = interior.inside_point();
// compute the intersection
halfspace_intersection_with_constructions_3(begin, end, P, origin, Traits());
// Compute the intersection of halfspaces by constructing explicitly
// the dual points with the default traits class for convex_hull_3.
// An interior point is not given.
template <class PlaneIterator, class Polyhedron>
void halfspace_intersection_with_constructions_3 (PlaneIterator pbegin,
PlaneIterator pend,
Polyhedron &P) {
halfspace_intersection_with_constructions_3(pbegin, pend, P, boost::none);
}
} // namespace CGAL