WFMath
1.0.2
|
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 AxisBox<dim> Ball<dim>::boundingBox() const 00040 { 00041 Point<dim> p_low, p_high; 00042 00043 for(int i = 0; i < dim; ++i) { 00044 p_low[i] = m_center[i] - m_radius; 00045 p_high[i] = m_center[i] + m_radius; 00046 } 00047 00048 bool valid = m_center.isValid(); 00049 00050 p_low.setValid(valid); 00051 p_high.setValid(valid); 00052 00053 return AxisBox<dim>(p_low, p_high, true); 00054 } 00055 00056 template<int dim, template<class, class> class container> 00057 Ball<dim> BoundingSphere(const container<Point<dim>, std::allocator<Point<dim> > >& c) 00058 { 00059 _miniball::Miniball<dim> m; 00060 _miniball::Wrapped_array<dim> w; 00061 00062 typename container<Point<dim>, std::allocator<Point<dim> > >::const_iterator i, end = c.end(); 00063 bool valid = true; 00064 00065 for(i = c.begin(); i != end; ++i) { 00066 valid = valid && i->isValid(); 00067 for(int j = 0; j < dim; ++j) 00068 w[j] = (*i)[j]; 00069 m.check_in(w); 00070 } 00071 00072 m.build(); 00073 00074 #ifndef NDEBUG 00075 double dummy; 00076 #endif 00077 assert("Check that bounding sphere is good to library accuracy" && 00078 m.accuracy(dummy) < numeric_constants<CoordType>::epsilon()); 00079 00080 w = m.center(); 00081 Point<dim> center; 00082 00083 for(int j = 0; j < dim; ++j) 00084 center[j] = w[j]; 00085 00086 center.setValid(valid); 00087 00088 return Ball<dim>(center, std::sqrt(m.squared_radius())); 00089 } 00090 00091 template<int dim, template<class, class> class container> 00092 Ball<dim> BoundingSphereSloppy(const container<Point<dim>, std::allocator<Point<dim> > >& c) 00093 { 00094 // This is based on the algorithm given by Jack Ritter 00095 // in Volume 2, Number 4 of Ray Tracing News 00096 // <http://www.acm.org/tog/resources/RTNews/html/rtnews7b.html> 00097 00098 typename container<Point<dim>, std::allocator<Point<dim> > >::const_iterator i = c.begin(), 00099 end = c.end(); 00100 if (i == end) { 00101 return Ball<dim>(); 00102 } 00103 00104 CoordType min[dim], max[dim]; 00105 typename container<Point<dim>, std::allocator<Point<dim> > >::const_iterator min_p[dim], max_p[dim]; 00106 bool valid = i->isValid(); 00107 00108 for(int j = 0; j < dim; ++j) { 00109 min[j] = max[j] = (*i)[j]; 00110 min_p[j] = max_p[j] = i; 00111 } 00112 00113 while(++i != end) { 00114 valid = valid && i->isValid(); 00115 for(int j = 0; j < dim; ++j) { 00116 if(min[j] > (*i)[j]) { 00117 min[j] = (*i)[j]; 00118 min_p[j] = i; 00119 } 00120 if(max[j] < (*i)[j]) { 00121 max[j] = (*i)[j]; 00122 max_p[j] = i; 00123 } 00124 } 00125 } 00126 00127 CoordType span = -1; 00128 int direction = -1; 00129 00130 for(int j = 0; j < dim; ++j) { 00131 CoordType new_span = max[j] - min[j]; 00132 if(new_span > span) { 00133 span = new_span; 00134 direction = j; 00135 } 00136 } 00137 00138 assert("Have a direction of maximum size" && direction != -1); 00139 00140 Point<dim> center = Midpoint(*(min_p[direction]), *(max_p[direction])); 00141 CoordType dist = SloppyDistance(*(min_p[direction]), center); 00142 00143 for(i = c.begin(); i != end; ++i) { 00144 if(i == min_p[direction] || i == max_p[direction]) 00145 continue; // We already have these 00146 00147 CoordType new_dist = SloppyDistance(*i, center); 00148 00149 if(new_dist > dist) { 00150 CoordType delta_dist = (new_dist - dist) / 2; 00151 // Even though new_dist may be too large, delta_dist / new_dist 00152 // always gives enough of a shift to include the new point. 00153 center += (*i - center) * delta_dist / new_dist; 00154 dist += delta_dist; 00155 assert("Shifted ball contains new point" && 00156 SquaredDistance(*i, center) <= dist * dist); 00157 } 00158 } 00159 00160 center.setValid(valid); 00161 00162 return Ball<dim>(center, dist); 00163 } 00164 00165 // These two are here, instead of defined in the class, to 00166 // avoid include order problems 00167 00168 template<int dim> 00169 inline Ball<dim> Point<dim>::boundingSphere() const 00170 { 00171 return Ball<dim>(*this, 0); 00172 } 00173 00174 template<int dim> 00175 inline Ball<dim> Point<dim>::boundingSphereSloppy() const 00176 { 00177 return Ball<dim>(*this, 0); 00178 } 00179 00180 } // namespace WFMath 00181 00182 #endif // WFMATH_BALL_FUNCS_H