curve.cpp Source File

00001 /*=============================================================================
00002         File: curve.cpp
00003      Purpose:       
00004     Revision: $Id: curve.cpp,v 1.3 2003/01/13 19:41:21 philosophil Exp $
00005   Created by: Philippe Lavoie          (3 Oct, 1996)
00006  Modified by: 
00007 
00008  Copyright notice:
00009           Copyright (C) 1996-1997 Philippe Lavoie
00010  
00011           This library is free software; you can redistribute it and/or
00012           modify it under the terms of the GNU Library General Public
00013           License as published by the Free Software Foundation; either
00014           version 2 of the License, or (at your option) any later version.
00015  
00016           This library is distributed in the hope that it will be useful,
00017           but WITHOUT ANY WARRANTY; without even the implied warranty of
00018           MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
00019           Library General Public License for more details.
00020  
00021           You should have received a copy of the GNU Library General Public
00022           License along with this library; if not, write to the Free
00023           Software Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA.
00024 =============================================================================*/
00025 
00026 #ifndef PLIB_NURBS_CURVE_SOURCE
00027 #define PLIB_NURBS_CURVE_SOURCE
00028 
00029 
00030 #include <curve.h>
00031 
00034 namespace PLib {
00035 
00063 template <class T, int N>
00064 T ParaCurve<T,N>::minDist2(const Point_nD<T,N>& p, T& guess,T error,T s,int sep,int maxIter, T um, T uM) const {
00065   if(um<0)
00066     um = minKnot() ;
00067   if(uM<0)
00068     uM = maxKnot() ;
00069   if(s<0)
00070     s = uM-um ;
00071   T d,d1,d2 ;
00072   Point_nD<T,N> p2 ;
00073   p2 = pointAt(guess) ;
00074   d = norm2(p-p2) ;
00075   d2 = d1 = 0 ;
00076   int niter = 0 ;
00077   T u1,u2 ;
00078   T step ;
00079   step = s/(T)sep ;
00080   u1 = guess-s ;
00081   u2 = guess+s ;
00082   while(d>error && niter<maxIter) {
00083     if(u1<um)
00084       u1=um;
00085     if(u2>uM)
00086       u2 = uM ;
00087     T u = u1 ;
00088     d2 = d1 ;
00089     for(;u<u2;u+=step){
00090       p2 = pointAt(u) ;
00091       d1 = norm2(p-p2) ;
00092       if(d1<d){
00093         d = d1 ;
00094         guess = u ;
00095       }
00096     }
00097     s /= 2.0 ;
00098     u1 = guess - s ;
00099     u2 = guess + s ;
00100     step = 2.0*s/(T)sep ;
00101     if(d-d2==0) niter = maxIter ;
00102     if(step<error) niter = maxIter ;
00103     niter++;
00104   }
00105   return d ;
00106 }
00107 
00108 
00137 template <class T, int N>
00138 Point_nD<T,N> ParaCurve<T,N>::minDistY(T y, T& guessU, T error, T s, int sep, int maxIter, T um, T uM) const {
00139   if(um<0)
00140     um = minKnot() ;
00141   if(uM<0)
00142     uM = maxKnot() ;
00143   if(s<0)
00144     s = uM-um ;
00145   T d,d1,d2 ;
00146   Point_nD<T,N> p2,result ;
00147   p2 = pointAt(guessU) ;
00148   result = p2 ;
00149   d = to2power(y-p2.y()) ;
00150   d2 = d1 = 0 ;
00151   int niter = 0 ;
00152   T u1,u2 ;
00153   T step ;
00154   step = s/(T)sep ;
00155   u1 = guessU-s ;
00156   u2 = guessU+s ;
00157   while(d>error && niter<maxIter) {
00158     if(u1<um)
00159       u1=um;
00160     if(u2>uM)
00161       u2 = uM ;
00162     T u = u1 ;
00163     d2 = d1 ;
00164     for(;u<u2;u+=step){
00165       p2 = pointAt(u) ;
00166       d1 = to2power(y-p2.y()) ;
00167       if(d1<d){
00168         d = d1 ;
00169         guessU = u ;
00170         result = p2 ;
00171       }
00172     }
00173     s /= 2.0 ;
00174     u1 = guessU - s ;
00175     u2 = guessU + s ;
00176     step = 2.0*s/(T)sep ;
00177     if(d-d2==0) niter = maxIter ;
00178     if(step<error) niter = maxIter ;
00179     niter++;
00180   }
00181   return result ;
00182 }
00183 
00184 
00213 template <class T, int N>
00214 Point_nD<T,N> ParaCurve<T,N>::minDistX(T x, T& guessU, T error, T s, int sep, int maxIter, T um, T uM) const {
00215   if(um<0)
00216     um = minKnot() ;
00217   if(uM<0)
00218     uM = maxKnot() ;
00219   if(s<0)
00220     s = uM-um ;
00221   T d,d1,d2 ;
00222   Point_nD<T,N> p2,result ;
00223   p2 = pointAt(guessU) ;
00224   result = p2 ;
00225   d = to2power(x-p2.x()) ;
00226   d2 = d1 = 0 ;
00227   int niter = 0 ;
00228   T u1,u2 ;
00229   T step ;
00230   step = s/(T)sep ;
00231   u1 = guessU-s ;
00232   u2 = guessU+s ;
00233   while(d>error && niter<maxIter) {
00234     if(u1<um)
00235       u1=um;
00236     if(u2>uM)
00237       u2 = uM ;
00238     T u = u1 ;
00239     d2 = d1 ;
00240     for(;u<u2;u+=step){
00241       p2 = pointAt(u) ;
00242       d1 = to2power(x-p2.x()) ;
00243       if(d1<d){
00244         d = d1 ;
00245         guessU = u ;
00246         result = p2 ;
00247       }
00248     }
00249     s /= 2.0 ;
00250     u1 = guessU - s ;
00251     u2 = guessU + s ;
00252     step = 2.0*s/(T)sep ;
00253     if(d-d2==0) niter = maxIter ;
00254     if(step<error) niter = maxIter ;
00255     niter++;
00256   }
00257   return result ;
00258 }
00259 
00260 
00289 template <class T, int N>
00290 Point_nD<T,N> ParaCurve<T,N>::minDistZ(T z, T& guessU, T error, T s, int sep, int maxIter, T um, T uM) const {
00291   if(um<0)
00292     um = minKnot() ;
00293   if(uM<0)
00294     uM = maxKnot() ;
00295   if(s<0)
00296     s = uM-um ;
00297   T d,d1,d2 ;
00298   Point_nD<T,N> p2,result ;
00299   p2 = pointAt(guessU) ;
00300   result = p2 ;
00301   d = to2power(z-p2.z()) ;
00302   d2 = d1 = 0 ;
00303   int niter = 0 ;
00304   T u1,u2 ;
00305   T step ;
00306   step = s/(T)sep ;
00307   u1 = guessU-s ;
00308   u2 = guessU+s ;
00309   while(d>error && niter<maxIter) {
00310     if(u1<um)
00311       u1=um;
00312     if(u2>uM)
00313       u2 = uM ;
00314     T u = u1 ;
00315     d2 = d1 ;
00316     for(;u<u2;u+=step){
00317       p2 = pointAt(u) ;
00318       d1 = to2power(z-p2.z()) ;
00319       if(d1<d){
00320         d = d1 ;
00321         guessU = u ;
00322         result = p2 ;
00323       }
00324     }
00325     s /= 2.0 ;
00326     u1 = guessU - s ;
00327     u2 = guessU + s ;
00328     step = 2.0*s/(T)sep ;
00329     if(d-d2==0) niter = maxIter ;
00330     if(step<error) niter = maxIter ;
00331     niter++;
00332   }
00333   return result ;
00334 }
00335 
00336 template <class T, int N>
00337 inline T coordValue(CoordinateType coord, const Point_nD<T,N>& p){
00338   switch(coord){
00339   case coordX: return p.x() ; break ;
00340   case coordY: return p.y() ; break ;
00341   case coordZ: return p.z() ; break ;
00342   }
00343   return 0.0 ; // elliminates warning messages
00344 }
00345 
00370 template <class T, int N>
00371 T ParaCurve<T,N>::extremum(int findMin, CoordinateType coord, T minDu, int sep, int maxIter, T um, T uM) const {
00372   if(um<0)
00373     um = minKnot() ;
00374   if(uM<0)
00375     uM = maxKnot() ;
00376 
00377   T c,d,d1,d2,result,guessU ;
00378   Point_nD<T,N> p2 ;
00379   p2 = pointAt(um) ;
00380   c = coordValue(coord,p2) ;
00381   p2 = pointAt(uM) ;
00382   if(findMin)
00383     c = minimum(c,coordValue(coord,p2)) ;
00384   else
00385     c = maximum(c,coordValue(coord,p2)) ;
00386   result = c ;
00387   d = minDu*10.0 ;
00388   d2 = d1 = 0 ;
00389   int niter = 0 ;
00390   T u1,u2 ;
00391   T step ;
00392   T s ;
00393   s = uM - um ;
00394   step = s/(T)(sep+1) ;
00395   u1 = um ;
00396   u2 = uM ;
00397   guessU = um ;
00398   while(d>minDu && niter<maxIter) {
00399     if(u1<um)
00400       u1=um;
00401     if(u2>uM)
00402       u2 = uM ;
00403     T u = u1 ;
00404     d2 = c ;
00405     d = guessU ;
00406     for(;u<=u2;u+=step){
00407       p2 = pointAt(u) ;
00408       if(findMin){
00409         d1 = minimum(c,coordValue(coord,p2)) ;
00410         if(d1<c){
00411           c = d1 ;
00412           guessU = u ;
00413           result = d1 ;
00414         }
00415       }
00416       else{
00417         d1 = maximum(c,coordValue(coord,p2)) ;
00418         if(d1>c){
00419           c = d1 ;
00420           guessU = u ;
00421           result = d1 ;
00422         }
00423       }
00424     }
00425     s /= 2.0 ;
00426     u1 = guessU - s ;
00427     u2 = guessU + s ;
00428     step = 2.0*s/(T)sep ;
00429     if((c-d2)==0.0) niter = maxIter ;
00430     if(step<minDu) niter = maxIter ;
00431     d = absolute(guessU-d) ;
00432     niter++;
00433   }
00434   return result ;
00435 }
00436 
00437 } // end namespace
00438 
00439 #endif // #define PLIB_NURBS_CURVE_SOURCE