rheolef  6.5
point_predicate.cc
Go to the documentation of this file.
1 // together witha custom cgal kernel that uses rheolef::point_basic<T>
2 // "CGAL/internal/Static_filters/Orientation_3.h"
3 #include "rheolef/point.h"
4 
5 #ifdef _RHEOLEF_HAVE_CGAL
6 #include "rheolef/cgal_traits.h"
7 #endif // _RHEOLEF_HAVE_CGAL
8 
9 namespace rheolef {
10 
11 template<class T>
12 static
13 T
15  const point_basic<T>& b)
16 {
17  T ax0 = a[0] - x[0];
18  T bx0 = b[0] - x[0];
19  T ax1 = a[1] - x[1];
20  T bx1 = b[1] - x[1];
21  return ax0*bx1 - ax1*bx0;
22 }
23 template <class T>
24 static
25 T
27  const point_basic<T>& b, const point_basic<T>& c)
28 {
29  T ax0 = a[0] - x[0];
30  T bx0 = b[0] - x[0];
31  T cx0 = c[0] - x[0];
32  T ax1 = a[1] - x[1];
33  T bx1 = b[1] - x[1];
34  T cx1 = c[1] - x[1];
35  T ax2 = a[2] - x[2];
36  T bx2 = b[2] - x[2];
37  T cx2 = c[2] - x[2];
38  return ax0 * (bx1 * cx2 - bx2 * cx1)
39  + bx0 * (cx1 * ax2 - cx2 * ax1)
40  + cx0 * (ax1 * bx2 - ax2 * bx1);
41 }
42 template <class T>
43 int
45  const point_basic<T>& a,
46  const point_basic<T>& b,
47  const point_basic<T>& c)
48 {
49 #ifdef _RHEOLEF_HAVE_CGAL
50  typedef typename geo_cgal_traits<T,2>::Kernel Kernel;
51  typename Kernel::Orientation_2 orientation;
52  CGAL::Orientation sgn = orientation(a, b, c);
53  return (sgn == CGAL::NEGATIVE) ? -1 : ((sgn == CGAL::ZERO) ? 0 : 1);
54 #else // _RHEOLEF_HAVE_CGAL
55  T sgn = inexact_orient2d(a, b, c);
56  return (sgn < 0) ? -1 : ((sgn == 0) ? 0 : 1);
57 #endif // _RHEOLEF_HAVE_CGAL
58 }
59 template <class T>
60 int
62  const point_basic<T>& a,
63  const point_basic<T>& b,
64  const point_basic<T>& c,
65  const point_basic<T>& d)
66 {
67 #ifdef _RHEOLEF_HAVE_CGAL
68  typedef typename geo_cgal_traits<T,3>::Kernel Kernel;
69  typename Kernel::Orientation_3 orientation;
70  CGAL::Orientation sgn = orientation(a, b, c, d);
71  return (sgn == CGAL::NEGATIVE) ? -1 : ((sgn == CGAL::ZERO) ? 0 : 1);
72 #else // _RHEOLEF_HAVE_CGAL
73  T sgn = inexact_orient3d(a, b, c, d);
74  return (sgn < 0) ? -1 : ((sgn == 0) ? 0 : 1);
75 #endif // _RHEOLEF_HAVE_CGAL
76 }
77 template<class T>
78 static
79 T
81  const point_basic<T>& c)
82 {
83  T value = inexact_orient2d(a, b, c);
84 #ifndef _RHEOLEF_HAVE_CGAL
85  return value;
86 #else // _RHEOLEF_HAVE_CGAL
87  int sgn = sign_orient2d(a, b, c);
88  value = fabs(value);
89  if (sgn == 0) return 0;
90  if (value != 0) return sgn*value;
91  // sgn != 0 but value == 0
93 #endif // _RHEOLEF_HAVE_CGAL
94 }
95 template <class T>
96 static
97 T
99  const point_basic<T>& c, const point_basic<T>& d)
100 {
101  T value = inexact_orient3d(a, b, c, d);
102 #ifndef _RHEOLEF_HAVE_CGAL
103  return value;
104 #else // _RHEOLEF_HAVE_CGAL
105  int sgn = sign_orient3d(a, b, c, d);
106  value = fabs(value);
107  if (sgn == 0) return 0;
108  if (value != 0) return sgn*value;
109  // sgn != 0 but value == 0
110  return sgn*std::numeric_limits<T>::epsilon();
111 #endif // _RHEOLEF_HAVE_CGAL
112 }
113 #define _RHEOLEF_instanciation(T) \
114 template T orient2d ( \
115  const point_basic<T>&, \
116  const point_basic<T>&, \
117  const point_basic<T>&); \
118 template T orient3d ( \
119  const point_basic<T>&, \
120  const point_basic<T>&, \
121  const point_basic<T>&, \
122  const point_basic<T>&); \
123 template int sign_orient2d ( \
124  const point_basic<T>&, \
125  const point_basic<T>&, \
126  const point_basic<T>&); \
127 template int sign_orient3d ( \
128  const point_basic<T>&, \
129  const point_basic<T>&, \
130  const point_basic<T>&, \
131  const point_basic<T>&); \
132 
134 
135 } // namespace rheolef
136