Loading...
Searching...
No Matches
point_predicate.cc
Go to the documentation of this file.
// robust floating point predicates:
// implementation using exact CGAL predicates, when available
// together witha custom cgal kernel that uses rheolef::point_basic<T>
//
// author: Pierre.Saramito@imag.fr
//
// date: 12 march 2012
//
// all the work for exact predicates is performed by :
// "CGAL/internal/Static_filters/Orientation_3.h"
//
#include "rheolef/point.h"
#ifdef _RHEOLEF_HAVE_CGAL
#include "rheolef/cgal_traits.h"
#endif // _RHEOLEF_HAVE_CGAL
namespace rheolef {
template<class T>
static
T
inexact_orient2d (const point_basic<T>& x, const point_basic<T>& a,
const point_basic<T>& b)
{
T ax0 = a[0] - x[0];
T bx0 = b[0] - x[0];
T ax1 = a[1] - x[1];
T bx1 = b[1] - x[1];
return ax0*bx1 - ax1*bx0;
}
template <class T>
static
T
inexact_orient3d (const point_basic<T>& x, const point_basic<T>& a,
const point_basic<T>& b, const point_basic<T>& c)
{
T ax0 = a[0] - x[0];
T bx0 = b[0] - x[0];
T cx0 = c[0] - x[0];
T ax1 = a[1] - x[1];
T bx1 = b[1] - x[1];
T cx1 = c[1] - x[1];
T ax2 = a[2] - x[2];
T bx2 = b[2] - x[2];
T cx2 = c[2] - x[2];
return ax0 * (bx1 * cx2 - bx2 * cx1)
+ bx0 * (cx1 * ax2 - cx2 * ax1)
+ cx0 * (ax1 * bx2 - ax2 * bx1);
}
template <class T>
int
sign_orient2d (
const point_basic<T>& a,
const point_basic<T>& b,
const point_basic<T>& c)
{
#ifdef _RHEOLEF_HAVE_CGAL
typedef typename geo_cgal_traits<T,2>::Kernel Kernel;
typename Kernel::Orientation_2 orientation;
CGAL::Orientation sgn = orientation(a, b, c);
return (sgn == CGAL::NEGATIVE) ? -1 : ((sgn == CGAL::ZERO) ? 0 : 1);
#else // _RHEOLEF_HAVE_CGAL
T sgn = inexact_orient2d(a, b, c);
return (sgn < 0) ? -1 : ((sgn == 0) ? 0 : 1);
#endif // _RHEOLEF_HAVE_CGAL
}
template <class T>
int
sign_orient3d (
const point_basic<T>& a,
const point_basic<T>& b,
const point_basic<T>& c,
const point_basic<T>& d)
{
#ifdef _RHEOLEF_HAVE_CGAL
typedef typename geo_cgal_traits<T,3>::Kernel Kernel;
typename Kernel::Orientation_3 orientation;
CGAL::Orientation sgn = orientation(a, b, c, d);
return (sgn == CGAL::NEGATIVE) ? -1 : ((sgn == CGAL::ZERO) ? 0 : 1);
#else // _RHEOLEF_HAVE_CGAL
T sgn = inexact_orient3d(a, b, c, d);
return (sgn < 0) ? -1 : ((sgn == 0) ? 0 : 1);
#endif // _RHEOLEF_HAVE_CGAL
}
template<class T>
T
orient2d (const point_basic<T>& a, const point_basic<T>& b,
const point_basic<T>& c)
{
T value = inexact_orient2d(a, b, c);
#ifndef _RHEOLEF_HAVE_CGAL
return value;
#else // _RHEOLEF_HAVE_CGAL
int sgn = sign_orient2d(a, b, c);
value = fabs(value);
if (sgn == 0) return 0;
if (value != 0) return sgn*value;
// sgn != 0 but value == 0
return sgn*std::numeric_limits<T>::epsilon();
#endif // _RHEOLEF_HAVE_CGAL
}
template <class T>
T
orient3d (const point_basic<T>& a, const point_basic<T>& b,
const point_basic<T>& c, const point_basic<T>& d)
{
T value = inexact_orient3d(a, b, c, d);
#ifndef _RHEOLEF_HAVE_CGAL
return value;
#else // _RHEOLEF_HAVE_CGAL
int sgn = sign_orient3d(a, b, c, d);
value = fabs(value);
if (sgn == 0) return 0;
if (value != 0) return sgn*value;
// sgn != 0 but value == 0
return sgn*std::numeric_limits<T>::epsilon();
#endif // _RHEOLEF_HAVE_CGAL
}
// ----------------------------------------------------------------------------
// instanciation in library
// ----------------------------------------------------------------------------
#define _RHEOLEF_instanciation(T) \
template T orient2d ( \
const point_basic<T>&, \
const point_basic<T>&, \
const point_basic<T>&); \
template T orient3d ( \
const point_basic<T>&, \
const point_basic<T>&, \
const point_basic<T>&, \
const point_basic<T>&); \
template int sign_orient2d ( \
const point_basic<T>&, \
const point_basic<T>&, \
const point_basic<T>&); \
template int sign_orient3d ( \
const point_basic<T>&, \
const point_basic<T>&, \
const point_basic<T>&, \
const point_basic<T>&); \
_RHEOLEF_instanciation(Float)
} // namespace rheolef
_RHEOLEF_instanciation
#define _RHEOLEF_instanciation(T, M, A)
Definition asr.cc:223
Float
see the Float page for the full documentation
rheolef::point_basic
Definition point.h:86
T
Expr1::float_type T
Definition field_expr.h:230
mkgeo_ball.b
int b
Definition mkgeo_ball.sh:152
mkgeo_ball.c
int c
Definition mkgeo_ball.sh:153
mkgeo_ball.a
int a
Definition mkgeo_ball.sh:151
rheolef
This file is part of Rheolef.
Definition compiler_eigen.h:39
rheolef::orient2d
T orient2d(const point_basic< T > &a, const point_basic< T > &b, const point_basic< T > &c)
Definition point_predicate.cc:107
rheolef::orient3d
T orient3d(const point_basic< T > &a, const point_basic< T > &b, const point_basic< T > &c, const point_basic< T > &d)
Definition point_predicate.cc:124
rheolef::sign_orient2d
int sign_orient2d(const point_basic< T > &a, const point_basic< T > &b, const point_basic< T > &c)
Definition point_predicate.cc:72
rheolef::sign_orient3d
int sign_orient3d(const point_basic< T > &a, const point_basic< T > &b, const point_basic< T > &c, const point_basic< T > &d)
Definition point_predicate.cc:89
d
size_t d
Definition neumann-laplace-lambda.cc:28
sgn
Float sgn(Float x)
Definition sgn.icc:25
rheolef::geo_cgal_traits
Definition cgal_traits.h:34