WFMath
0.3.12
|
00001 // ball_funcs.h (n-dimensional ball implementation) 00002 // 00003 // The WorldForge Project 00004 // Copyright (C) 2000, 2001 The WorldForge Project 00005 // 00006 // This program is free software; you can redistribute it and/or modify 00007 // it under the terms of the GNU General Public License as published by 00008 // the Free Software Foundation; either version 2 of the License, or 00009 // (at your option) any later version. 00010 // 00011 // This program is distributed in the hope that it will be useful, 00012 // but WITHOUT ANY WARRANTY; without even the implied warranty of 00013 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 00014 // GNU General Public License for more details. 00015 // 00016 // You should have received a copy of the GNU General Public License 00017 // along with this program; if not, write to the Free Software 00018 // Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA. 00019 // 00020 // For information about WorldForge and its authors, please contact 00021 // the Worldforge Web Site at http://www.worldforge.org. 00022 // 00023 00024 // Author: Ron Steinke 00025 00026 #ifndef WFMATH_BALL_FUNCS_H 00027 #define WFMATH_BALL_FUNCS_H 00028 00029 #include <wfmath/ball.h> 00030 00031 #include <wfmath/axisbox.h> 00032 #include <wfmath/miniball.h> 00033 00034 #include <cassert> 00035 00036 namespace WFMath { 00037 00038 template<int dim> 00039 inline bool Ball<dim>::isEqualTo(const Ball<dim>& b, double epsilon) const 00040 { 00041 return Equal(m_center, b.m_center, epsilon) 00042 && Equal(m_radius, b.m_radius, epsilon); 00043 } 00044 00045 template<int dim> 00046 AxisBox<dim> Ball<dim>::boundingBox() const 00047 { 00048 Point<dim> p_low, p_high; 00049 00050 for(int i = 0; i < dim; ++i) { 00051 p_low[i] = m_center[i] - m_radius; 00052 p_high[i] = m_center[i] + m_radius; 00053 } 00054 00055 bool valid = m_center.isValid(); 00056 00057 p_low.setValid(valid); 00058 p_high.setValid(valid); 00059 00060 return AxisBox<dim>(p_low, p_high, true); 00061 } 00062 00063 template<int dim, template<class, class> class container> 00064 Ball<dim> BoundingSphere(const container<Point<dim>, std::allocator<Point<dim> > >& c) 00065 { 00066 _miniball::Miniball<dim> m; 00067 _miniball::Wrapped_array<dim> w; 00068 00069 typename container<Point<dim>, std::allocator<Point<dim> > >::const_iterator i, end = c.end(); 00070 bool valid = true; 00071 00072 for(i = c.begin(); i != end; ++i) { 00073 valid = valid && i->isValid(); 00074 for(int j = 0; j < dim; ++j) 00075 w[j] = (*i)[j]; 00076 m.check_in(w); 00077 } 00078 00079 m.build(); 00080 00081 #ifndef NDEBUG 00082 double dummy; 00083 #endif 00084 assert("Check that bounding sphere is good to library accuracy" && 00085 m.accuracy(dummy) < WFMATH_EPSILON); 00086 00087 w = m.center(); 00088 Point<dim> center; 00089 00090 for(int j = 0; j < dim; ++j) 00091 center[j] = w[j]; 00092 00093 center.setValid(valid); 00094 00095 return Ball<dim>(center, std::sqrt(m.squared_radius())); 00096 } 00097 00098 template<int dim, template<class, class> class container> 00099 Ball<dim> BoundingSphereSloppy(const container<Point<dim>, std::allocator<Point<dim> > >& c) 00100 { 00101 // This is based on the algorithm given by Jack Ritter 00102 // in Volume 2, Number 4 of Ray Tracing News 00103 // <http://www.acm.org/tog/resources/RTNews/html/rtnews7b.html> 00104 00105 typename container<Point<dim>, std::allocator<Point<dim> > >::const_iterator i = c.begin(), 00106 end = c.end(); 00107 if (i == end) { 00108 return Ball<dim>(); 00109 } 00110 00111 CoordType min[dim], max[dim]; 00112 typename container<Point<dim>, std::allocator<Point<dim> > >::const_iterator min_p[dim], max_p[dim]; 00113 bool valid = i->isValid(); 00114 00115 for(int j = 0; j < dim; ++j) { 00116 min[j] = max[j] = (*i)[j]; 00117 min_p[j] = max_p[j] = i; 00118 } 00119 00120 while(++i != end) { 00121 valid = valid && i->isValid(); 00122 for(int j = 0; j < dim; ++j) { 00123 if(min[j] > (*i)[j]) { 00124 min[j] = (*i)[j]; 00125 min_p[j] = i; 00126 } 00127 if(max[j] < (*i)[j]) { 00128 max[j] = (*i)[j]; 00129 max_p[j] = i; 00130 } 00131 } 00132 } 00133 00134 CoordType span = -1; 00135 int direction = -1; 00136 00137 for(int j = 0; j < dim; ++j) { 00138 CoordType new_span = max[j] - min[j]; 00139 if(new_span > span) { 00140 span = new_span; 00141 direction = j; 00142 } 00143 } 00144 00145 assert("Have a direction of maximum size" && direction != -1); 00146 00147 Point<dim> center = Midpoint(*(min_p[direction]), *(max_p[direction])); 00148 CoordType dist = SloppyDistance(*(min_p[direction]), center); 00149 00150 for(i = c.begin(); i != end; ++i) { 00151 if(i == min_p[direction] || i == max_p[direction]) 00152 continue; // We already have these 00153 00154 CoordType new_dist = SloppyDistance(*i, center); 00155 00156 if(new_dist > dist) { 00157 CoordType delta_dist = (new_dist - dist) / 2; 00158 // Even though new_dist may be too large, delta_dist / new_dist 00159 // always gives enough of a shift to include the new point. 00160 center += (*i - center) * delta_dist / new_dist; 00161 dist += delta_dist; 00162 assert("Shifted ball contains new point" && 00163 SquaredDistance(*i, center) <= dist * dist); 00164 } 00165 } 00166 00167 center.setValid(valid); 00168 00169 return Ball<2>(center, dist); 00170 } 00171 00172 // These two are here, instead of defined in the class, to 00173 // avoid include order problems 00174 00175 template<int dim> 00176 inline Ball<dim> Point<dim>::boundingSphere() const 00177 { 00178 return Ball<dim>(*this, 0); 00179 } 00180 00181 template<int dim> 00182 inline Ball<dim> Point<dim>::boundingSphereSloppy() const 00183 { 00184 return Ball<dim>(*this, 0); 00185 } 00186 00187 } // namespace WFMath 00188 00189 #endif // WFMATH_BALL_FUNCS_H