31#include "rheolef/point.h"
33#ifdef _RHEOLEF_HAVE_CGAL
34#include "rheolef/cgal_traits.h"
42inexact_orient2d (
const point_basic<T>& x,
const point_basic<T>& a,
43 const point_basic<T>& b)
49 return ax0*bx1 - ax1*bx0;
54inexact_orient3d (
const point_basic<T>& x,
const point_basic<T>& a,
55 const point_basic<T>& b,
const point_basic<T>& c)
66 return ax0 * (bx1 * cx2 - bx2 * cx1)
67 + bx0 * (cx1 * ax2 - cx2 * ax1)
68 + cx0 * (ax1 * bx2 - ax2 * bx1);
77#ifdef _RHEOLEF_HAVE_CGAL
79 typename Kernel::Orientation_2 orientation;
80 CGAL::Orientation
sgn = orientation(a, b, c);
81 return (
sgn == CGAL::NEGATIVE) ? -1 : ((
sgn == CGAL::ZERO) ? 0 : 1);
83 T sgn = inexact_orient2d(a, b, c);
84 return (
sgn < 0) ? -1 : ((
sgn == 0) ? 0 : 1);
95#ifdef _RHEOLEF_HAVE_CGAL
97 typename Kernel::Orientation_3 orientation;
98 CGAL::Orientation
sgn = orientation(a, b, c,
d);
99 return (
sgn == CGAL::NEGATIVE) ? -1 : ((
sgn == CGAL::ZERO) ? 0 : 1);
101 T sgn = inexact_orient3d(a, b, c,
d);
102 return (
sgn < 0) ? -1 : ((
sgn == 0) ? 0 : 1);
110 T value = inexact_orient2d(a, b, c);
111#ifndef _RHEOLEF_HAVE_CGAL
116 if (
sgn == 0)
return 0;
117 if (value != 0)
return sgn*value;
119 return sgn*std::numeric_limits<T>::epsilon();
127 T value = inexact_orient3d(a, b, c,
d);
128#ifndef _RHEOLEF_HAVE_CGAL
133 if (
sgn == 0)
return 0;
134 if (value != 0)
return sgn*value;
136 return sgn*std::numeric_limits<T>::epsilon();
142#define _RHEOLEF_instanciation(T) \
143template T orient2d ( \
144 const point_basic<T>&, \
145 const point_basic<T>&, \
146 const point_basic<T>&); \
147template T orient3d ( \
148 const point_basic<T>&, \
149 const point_basic<T>&, \
150 const point_basic<T>&, \
151 const point_basic<T>&); \
152template int sign_orient2d ( \
153 const point_basic<T>&, \
154 const point_basic<T>&, \
155 const point_basic<T>&); \
156template int sign_orient3d ( \
157 const point_basic<T>&, \
158 const point_basic<T>&, \
159 const point_basic<T>&, \
160 const point_basic<T>&); \
#define _RHEOLEF_instanciation(T, M, A)
see the Float page for the full documentation
This file is part of Rheolef.
T orient2d(const point_basic< T > &a, const point_basic< T > &b, const point_basic< T > &c)
T orient3d(const point_basic< T > &a, const point_basic< T > &b, const point_basic< T > &c, const point_basic< T > &d)
int sign_orient2d(const point_basic< T > &a, const point_basic< T > &b, const point_basic< T > &c)
int sign_orient3d(const point_basic< T > &a, const point_basic< T > &b, const point_basic< T > &c, const point_basic< T > &d)