WFMath  0.3.12
ball_funcs.h
1 // ball_funcs.h (n-dimensional ball implementation)
2 //
3 // The WorldForge Project
4 // Copyright (C) 2000, 2001 The WorldForge Project
5 //
6 // This program is free software; you can redistribute it and/or modify
7 // it under the terms of the GNU General Public License as published by
8 // the Free Software Foundation; either version 2 of the License, or
9 // (at your option) any later version.
10 //
11 // This program is distributed in the hope that it will be useful,
12 // but WITHOUT ANY WARRANTY; without even the implied warranty of
13 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14 // GNU General Public License for more details.
15 //
16 // You should have received a copy of the GNU General Public License
17 // along with this program; if not, write to the Free Software
18 // Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA.
19 //
20 // For information about WorldForge and its authors, please contact
21 // the Worldforge Web Site at http://www.worldforge.org.
22 //
23 
24 // Author: Ron Steinke
25 
26 #ifndef WFMATH_BALL_FUNCS_H
27 #define WFMATH_BALL_FUNCS_H
28 
29 #include <wfmath/ball.h>
30 
31 #include <wfmath/axisbox.h>
32 #include <wfmath/miniball.h>
33 
34 #include <cassert>
35 
36 namespace WFMath {
37 
38 template<int dim>
39 inline bool Ball<dim>::isEqualTo(const Ball<dim>& b, double epsilon) const
40 {
41  return Equal(m_center, b.m_center, epsilon)
42  && Equal(m_radius, b.m_radius, epsilon);
43 }
44 
45 template<int dim>
46 AxisBox<dim> Ball<dim>::boundingBox() const
47 {
48  Point<dim> p_low, p_high;
49 
50  for(int i = 0; i < dim; ++i) {
51  p_low[i] = m_center[i] - m_radius;
52  p_high[i] = m_center[i] + m_radius;
53  }
54 
55  bool valid = m_center.isValid();
56 
57  p_low.setValid(valid);
58  p_high.setValid(valid);
59 
60  return AxisBox<dim>(p_low, p_high, true);
61 }
62 
63 template<int dim, template<class, class> class container>
64 Ball<dim> BoundingSphere(const container<Point<dim>, std::allocator<Point<dim> > >& c)
65 {
66  _miniball::Miniball<dim> m;
67  _miniball::Wrapped_array<dim> w;
68 
69  typename container<Point<dim>, std::allocator<Point<dim> > >::const_iterator i, end = c.end();
70  bool valid = true;
71 
72  for(i = c.begin(); i != end; ++i) {
73  valid = valid && i->isValid();
74  for(int j = 0; j < dim; ++j)
75  w[j] = (*i)[j];
76  m.check_in(w);
77  }
78 
79  m.build();
80 
81 #ifndef NDEBUG
82  double dummy;
83 #endif
84  assert("Check that bounding sphere is good to library accuracy" &&
85  m.accuracy(dummy) < WFMATH_EPSILON);
86 
87  w = m.center();
88  Point<dim> center;
89 
90  for(int j = 0; j < dim; ++j)
91  center[j] = w[j];
92 
93  center.setValid(valid);
94 
95  return Ball<dim>(center, std::sqrt(m.squared_radius()));
96 }
97 
98 template<int dim, template<class, class> class container>
99 Ball<dim> BoundingSphereSloppy(const container<Point<dim>, std::allocator<Point<dim> > >& c)
100 {
101  // This is based on the algorithm given by Jack Ritter
102  // in Volume 2, Number 4 of Ray Tracing News
103  // <http://www.acm.org/tog/resources/RTNews/html/rtnews7b.html>
104 
105  typename container<Point<dim>, std::allocator<Point<dim> > >::const_iterator i = c.begin(),
106  end = c.end();
107  if (i == end) {
108  return Ball<dim>();
109  }
110 
111  CoordType min[dim], max[dim];
112  typename container<Point<dim>, std::allocator<Point<dim> > >::const_iterator min_p[dim], max_p[dim];
113  bool valid = i->isValid();
114 
115  for(int j = 0; j < dim; ++j) {
116  min[j] = max[j] = (*i)[j];
117  min_p[j] = max_p[j] = i;
118  }
119 
120  while(++i != end) {
121  valid = valid && i->isValid();
122  for(int j = 0; j < dim; ++j) {
123  if(min[j] > (*i)[j]) {
124  min[j] = (*i)[j];
125  min_p[j] = i;
126  }
127  if(max[j] < (*i)[j]) {
128  max[j] = (*i)[j];
129  max_p[j] = i;
130  }
131  }
132  }
133 
134  CoordType span = -1;
135  int direction = -1;
136 
137  for(int j = 0; j < dim; ++j) {
138  CoordType new_span = max[j] - min[j];
139  if(new_span > span) {
140  span = new_span;
141  direction = j;
142  }
143  }
144 
145  assert("Have a direction of maximum size" && direction != -1);
146 
147  Point<dim> center = Midpoint(*(min_p[direction]), *(max_p[direction]));
148  CoordType dist = SloppyDistance(*(min_p[direction]), center);
149 
150  for(i = c.begin(); i != end; ++i) {
151  if(i == min_p[direction] || i == max_p[direction])
152  continue; // We already have these
153 
154  CoordType new_dist = SloppyDistance(*i, center);
155 
156  if(new_dist > dist) {
157  CoordType delta_dist = (new_dist - dist) / 2;
158  // Even though new_dist may be too large, delta_dist / new_dist
159  // always gives enough of a shift to include the new point.
160  center += (*i - center) * delta_dist / new_dist;
161  dist += delta_dist;
162  assert("Shifted ball contains new point" &&
163  SquaredDistance(*i, center) <= dist * dist);
164  }
165  }
166 
167  center.setValid(valid);
168 
169  return Ball<2>(center, dist);
170 }
171 
172 // These two are here, instead of defined in the class, to
173 // avoid include order problems
174 
175 template<int dim>
176 inline Ball<dim> Point<dim>::boundingSphere() const
177 {
178  return Ball<dim>(*this, 0);
179 }
180 
181 template<int dim>
182 inline Ball<dim> Point<dim>::boundingSphereSloppy() const
183 {
184  return Ball<dim>(*this, 0);
185 }
186 
187 } // namespace WFMath
188 
189 #endif // WFMATH_BALL_FUNCS_H
bool Equal(const C &c1, const C &c2, double epsilon=WFMATH_EPSILON)
Test for equality up to precision epsilon.
Definition: const.h:103
Ball< dim > BoundingSphereSloppy(const container< Point< dim >, std::allocator< Point< dim > > > &c)
get a bounding sphere for a set of points
Definition: ball_funcs.h:99
void setValid(bool valid=true)
make isValid() return true if you've initialized the point by hand
Definition: point.h:129
Point< dim > Midpoint(const Point< dim > &p1, const Point< dim > &p2, CoordType dist=0.5)
Definition: point_funcs.h:270
float CoordType
Basic floating point type.
Definition: const.h:79
Ball< dim > BoundingSphere(const container< Point< dim >, std::allocator< Point< dim > > > &c)
get the minimal bounding sphere for a set of points
Definition: ball_funcs.h:64
A dim dimensional point.
Definition: const.h:50
A dim dimensional ball.
Definition: ball.h:34