0029769: Uninitialized data with BSplCLib_Cache, BSplSLib_Cache
authorabv <abv@opencascade.com>
Sun, 10 Jun 2018 19:40:12 +0000 (22:40 +0300)
committerbugmaster <bugmaster@opencascade.com>
Fri, 20 Jul 2018 13:59:22 +0000 (16:59 +0300)
Implementation of classes BSplCLib_Cache and BSplSLib_Cache is revised:
- Common functionality dealing with spans along one parametric direction is separated to new struct BSplCLib_CacheParams
- Empty constructors are removed; copying is prohibited
- Code reconsidering degree and other parameters on each call to BuildCache() is eliminated; curve parameters must be the same in constructor and all calls to BuildCache()
- Extra call to BuildCache() from constructor is eliminated

src/BSplCLib/BSplCLib_Cache.cxx
src/BSplCLib/BSplCLib_Cache.hxx
src/BSplCLib/BSplCLib_CacheParams.hxx [new file with mode: 0644]
src/BSplCLib/FILES
src/BSplSLib/BSplSLib_Cache.cxx
src/BSplSLib/BSplSLib_Cache.hxx
src/Geom2dAdaptor/Geom2dAdaptor_Curve.cxx
src/GeomAdaptor/GeomAdaptor_Curve.cxx
src/GeomAdaptor/GeomAdaptor_Surface.cxx

index 59efb2f..f835c9b 100644 (file)
@@ -31,159 +31,71 @@ static Standard_Real* ConvertArray(const Handle(TColStd_HArray2OfReal)& theHArra
   return (Standard_Real*) &(anArray(anArray.LowerRow(), anArray.LowerCol()));
 }
 
-
-BSplCLib_Cache::BSplCLib_Cache()
-{
-  myPolesWeights.Nullify();
-  myIsRational = Standard_False;
-  mySpanStart = 0.0;
-  mySpanLength = 0.0;
-  mySpanIndex = 0;
-  myDegree = 0;
-  myFlatKnots.Nullify();
-}
-
 BSplCLib_Cache::BSplCLib_Cache(const Standard_Integer&        theDegree,
                                const Standard_Boolean&        thePeriodic,
                                const TColStd_Array1OfReal&    theFlatKnots,
-                               const TColgp_Array1OfPnt2d&    thePoles2d,
+                               const TColgp_Array1OfPnt2d&    /* only used to distinguish from 3d variant */,
                                const TColStd_Array1OfReal*    theWeights)
+: myIsRational(theWeights != NULL),
+  myParams (theDegree, thePeriodic, theFlatKnots)
 {
-  Standard_Real aCacheParam = theFlatKnots.Value(theFlatKnots.Lower() + theDegree);
-  BuildCache(aCacheParam, theDegree, thePeriodic, 
-             theFlatKnots, thePoles2d, theWeights);
+  Standard_Integer aPWColNumber = (myIsRational ? 3 : 2);
+  myPolesWeights = new TColStd_HArray2OfReal (1, theDegree + 1, 1, aPWColNumber);
 }
 
 BSplCLib_Cache::BSplCLib_Cache(const Standard_Integer&        theDegree,
                                const Standard_Boolean&        thePeriodic,
                                const TColStd_Array1OfReal&    theFlatKnots,
-                               const TColgp_Array1OfPnt&      thePoles,
+                               const TColgp_Array1OfPnt&      /* only used to distinguish from 2d variant */,
                                const TColStd_Array1OfReal*    theWeights)
+: myIsRational(theWeights != NULL),
+  myParams (theDegree, thePeriodic, theFlatKnots)
 {
-  Standard_Real aCacheParam = theFlatKnots.Value(theFlatKnots.Lower() + theDegree);
-  BuildCache(aCacheParam, theDegree, thePeriodic, 
-             theFlatKnots, thePoles, theWeights);
+  Standard_Integer aPWColNumber = (myIsRational ? 4 : 3);
+  myPolesWeights = new TColStd_HArray2OfReal (1, theDegree + 1, 1, aPWColNumber);
 }
 
-
 Standard_Boolean BSplCLib_Cache::IsCacheValid(Standard_Real theParameter) const
 {
-  Standard_Real aNewParam = theParameter;
-  if (!myFlatKnots.IsNull())
-    PeriodicNormalization(myFlatKnots->Array1(), aNewParam);
-
-  Standard_Real aDelta = aNewParam - mySpanStart;
-  return ((aDelta >= 0.0 || mySpanIndex == mySpanIndexMin) &&
-          (aDelta < mySpanLength || mySpanIndex == mySpanIndexMax));
+  return myParams.IsCacheValid (theParameter);
 }
 
-void BSplCLib_Cache::PeriodicNormalization(const TColStd_Array1OfReal& theFlatKnots, 
-                                           Standard_Real& theParameter) const
-{
-  Standard_Real aPeriod = theFlatKnots.Value(theFlatKnots.Upper() - myDegree) - 
-                          theFlatKnots.Value(myDegree + 1) ;
-  if (theParameter < theFlatKnots.Value(myDegree + 1))
-  {
-    Standard_Real aScale = IntegerPart(
-        (theFlatKnots.Value(myDegree + 1) - theParameter) / aPeriod);
-    theParameter += aPeriod * (aScale + 1.0);
-  }
-  if (theParameter > theFlatKnots.Value(theFlatKnots.Upper() - myDegree))
-  {
-    Standard_Real aScale = IntegerPart(
-        (theParameter - theFlatKnots.Value(theFlatKnots.Upper() - myDegree)) / aPeriod);
-    theParameter -= aPeriod * (aScale + 1.0);
-  }
-}
-
-
 void BSplCLib_Cache::BuildCache(const Standard_Real&           theParameter,
-                                const Standard_Integer&        theDegree,
-                                const Standard_Boolean&        thePeriodic,
                                 const TColStd_Array1OfReal&    theFlatKnots,
                                 const TColgp_Array1OfPnt2d&    thePoles2d,
                                 const TColStd_Array1OfReal*    theWeights)
 {
   // Normalize theParameter for periodical B-splines
-  Standard_Real aNewParam = theParameter;
-  if (thePeriodic)
-  {
-    PeriodicNormalization(theFlatKnots, aNewParam);
-    myFlatKnots = new TColStd_HArray1OfReal(1, theFlatKnots.Length());
-    myFlatKnots->ChangeArray1() = theFlatKnots;
-  }
-  else if (!myFlatKnots.IsNull()) // Periodical curve became non-periodical
-    myFlatKnots.Nullify();
-
-  // Change the size of cached data if needed
-  myIsRational = (theWeights != NULL);
-  Standard_Integer aPWColNumber = myIsRational ? 3 : 2;
-  if (theDegree > myDegree)
-    myPolesWeights = new TColStd_HArray2OfReal(1, theDegree + 1, 1, aPWColNumber);
-
-  myDegree = theDegree;
-  mySpanIndex = 0;
-  BSplCLib::LocateParameter(theDegree, theFlatKnots, BSplCLib::NoMults(), 
-                            aNewParam, thePeriodic, mySpanIndex, aNewParam);
-  mySpanStart  = theFlatKnots.Value(mySpanIndex);
-  mySpanLength = theFlatKnots.Value(mySpanIndex + 1) - mySpanStart;
-  mySpanIndexMin = thePeriodic ? 0 : myDegree + 1;
-  mySpanIndexMax = theFlatKnots.Length() - 1 - theDegree;
+  Standard_Real aNewParam = myParams.PeriodicNormalization (theParameter);
+  myParams.LocateParameter (aNewParam, theFlatKnots);
 
   // Calculate new cache data
-  BSplCLib::BuildCache(mySpanStart, mySpanLength, thePeriodic, theDegree, 
-                       mySpanIndex, theFlatKnots, thePoles2d, theWeights, 
-                       myPolesWeights->ChangeArray2());
+  BSplCLib::BuildCache (myParams.SpanStart, myParams.SpanLength, myParams.IsPeriodic,
+                        myParams.Degree, myParams.SpanIndex, theFlatKnots, thePoles2d,
+                        theWeights, myPolesWeights->ChangeArray2());
 }
 
 void BSplCLib_Cache::BuildCache(const Standard_Real&           theParameter,
-                                const Standard_Integer&        theDegree,
-                                const Standard_Boolean&        thePeriodic,
                                 const TColStd_Array1OfReal&    theFlatKnots,
                                 const TColgp_Array1OfPnt&      thePoles,
                                 const TColStd_Array1OfReal*    theWeights)
 {
   // Create list of knots with repetitions and normalize theParameter for periodical B-splines
-  Standard_Real aNewParam = theParameter;
-  if (thePeriodic)
-  {
-    PeriodicNormalization(theFlatKnots, aNewParam);
-    myFlatKnots = new TColStd_HArray1OfReal(1, theFlatKnots.Length());
-    myFlatKnots->ChangeArray1() = theFlatKnots;
-  }
-  else if (!myFlatKnots.IsNull()) // Periodical curve became non-periodical
-    myFlatKnots.Nullify();
-
-  // Change the size of cached data if needed
-  myIsRational = (theWeights != NULL);
-  Standard_Integer aPWColNumber = myIsRational ? 4 : 3;
-  if (theDegree > myDegree)
-    myPolesWeights = new TColStd_HArray2OfReal(1, theDegree + 1, 1, aPWColNumber);
-
-  myDegree = theDegree;
-  mySpanIndex = 0;
-  BSplCLib::LocateParameter(theDegree, theFlatKnots, BSplCLib::NoMults(), 
-                            aNewParam, thePeriodic, mySpanIndex, aNewParam);
-  mySpanStart  = theFlatKnots.Value(mySpanIndex);
-  mySpanLength = theFlatKnots.Value(mySpanIndex + 1) - mySpanStart;
-  mySpanIndexMin = thePeriodic ? 0 : myDegree + 1;
-  mySpanIndexMax = theFlatKnots.Length() - 1 - theDegree;
+  Standard_Real aNewParam = myParams.PeriodicNormalization (theParameter);
+  myParams.LocateParameter (aNewParam, theFlatKnots);
 
   // Calculate new cache data
-  BSplCLib::BuildCache(mySpanStart, mySpanLength, thePeriodic, theDegree
-                       mySpanIndex, theFlatKnots, thePoles, theWeights, 
-                       myPolesWeights->ChangeArray2());
+  BSplCLib::BuildCache (myParams.SpanStart, myParams.SpanLength, myParams.IsPeriodic
+                        myParams.Degree, myParams.SpanIndex, theFlatKnots, thePoles,
+                        theWeights, myPolesWeights->ChangeArray2());
 }
 
-
 void BSplCLib_Cache::CalculateDerivative(const Standard_Real&    theParameter, 
                                          const Standard_Integer& theDerivative, 
                                                Standard_Real&    theDerivArray) const
 {
-  Standard_Real aNewParameter = theParameter;
-  if (!myFlatKnots.IsNull()) // B-spline is periodical
-    PeriodicNormalization(myFlatKnots->Array1(), aNewParameter);
-  aNewParameter = (aNewParameter - mySpanStart) / mySpanLength;
+  Standard_Real aNewParameter = myParams.PeriodicNormalization (theParameter);
+  aNewParameter = (aNewParameter - myParams.SpanStart) / myParams.SpanLength;
 
   Standard_Real* aPolesArray = ConvertArray(myPolesWeights);
   Standard_Integer aDimension = myPolesWeights->RowLength(); // number of columns
@@ -199,23 +111,23 @@ void BSplCLib_Cache::CalculateDerivative(const Standard_Real&    theParameter,
   // When the degree of curve is lesser than the requested derivative,
   // nullify array cells corresponding to greater derivatives
   Standard_Integer aDerivative = theDerivative;
-  if (myDegree < theDerivative)
+  if (myParams.Degree < theDerivative)
   {
-    aDerivative = myDegree;
-    for (Standard_Integer ind = myDegree * aDimension; ind < (theDerivative + 1) * aDimension; ind++)
+    aDerivative = myParams.Degree;
+    for (Standard_Integer ind = myParams.Degree * aDimension; ind < (theDerivative + 1) * aDimension; ind++)
     {
       aPntDeriv[ind] = 0.0;
       (&theDerivArray)[ind] = 0.0; // should be cleared separately, because aPntDeriv may look to another memory area
     }
   }
 
-  PLib::EvalPolynomial(aNewParameter, aDerivative, myDegree, aDimension, 
+  PLib::EvalPolynomial(aNewParameter, aDerivative, myParams.Degree, aDimension, 
                        aPolesArray[0], aPntDeriv[0]);
   // Unnormalize derivatives since those are computed normalized
   Standard_Real aFactor = 1.0;
   for (Standard_Integer deriv = 1; deriv <= aDerivative; deriv++)
   {
-    aFactor /= mySpanLength;
+    aFactor /= myParams.SpanLength;
     for (Standard_Integer ind = 0; ind < aDimension; ind++)
       aPntDeriv[aDimension * deriv + ind] *= aFactor;
   }
@@ -227,17 +139,15 @@ void BSplCLib_Cache::CalculateDerivative(const Standard_Real&    theParameter,
 
 void BSplCLib_Cache::D0(const Standard_Real& theParameter, gp_Pnt2d& thePoint) const
 {
-  Standard_Real aNewParameter = theParameter;
-  if (!myFlatKnots.IsNull()) // B-spline is periodical
-    PeriodicNormalization(myFlatKnots->Array1(), aNewParameter);
-  aNewParameter = (aNewParameter - mySpanStart) / mySpanLength;
+  Standard_Real aNewParameter = myParams.PeriodicNormalization (theParameter);
+  aNewParameter = (aNewParameter - myParams.SpanStart) / myParams.SpanLength;
 
   Standard_Real* aPolesArray = ConvertArray(myPolesWeights);
   Standard_Real aPoint[4];
   Standard_Integer aDimension = myPolesWeights->RowLength(); // number of columns
 
-  PLib::NoDerivativeEvalPolynomial(aNewParameter, myDegree,
-                                   aDimension, myDegree * aDimension,
+  PLib::NoDerivativeEvalPolynomial(aNewParameter, myParams.Degree,
+                                   aDimension, myParams.Degree * aDimension,
                                    aPolesArray[0], aPoint[0]);
 
   thePoint.SetCoord(aPoint[0], aPoint[1]);
@@ -247,17 +157,15 @@ void BSplCLib_Cache::D0(const Standard_Real& theParameter, gp_Pnt2d& thePoint) c
 
 void BSplCLib_Cache::D0(const Standard_Real& theParameter, gp_Pnt& thePoint) const
 {
-  Standard_Real aNewParameter = theParameter;
-  if (!myFlatKnots.IsNull()) // B-spline is periodical
-    PeriodicNormalization(myFlatKnots->Array1(), aNewParameter);
-  aNewParameter = (aNewParameter - mySpanStart) / mySpanLength;
+  Standard_Real aNewParameter = myParams.PeriodicNormalization (theParameter);
+  aNewParameter = (aNewParameter - myParams.SpanStart) / myParams.SpanLength;
 
   Standard_Real* aPolesArray = ConvertArray(myPolesWeights);
   Standard_Real aPoint[4];
   Standard_Integer aDimension = myPolesWeights->RowLength(); // number of columns
 
-  PLib::NoDerivativeEvalPolynomial(aNewParameter, myDegree,
-                                   aDimension, myDegree * aDimension,
+  PLib::NoDerivativeEvalPolynomial(aNewParameter, myParams.Degree,
+                                   aDimension, myParams.Degree * aDimension,
                                    aPolesArray[0], aPoint[0]);
 
   thePoint.SetCoord(aPoint[0], aPoint[1], aPoint[2]);
index 74bcd14..a4bfd46 100644 (file)
@@ -19,7 +19,6 @@
 #include <Standard_Type.hxx>
 #include <Standard_Transient.hxx>
 
-
 #include <gp_Pnt2d.hxx>
 #include <gp_Pnt.hxx>
 #include <gp_Vec2d.hxx>
@@ -31,6 +30,8 @@
 #include <TColgp_Array1OfPnt.hxx>
 #include <TColgp_Array1OfPnt2d.hxx>
 
+#include <BSplCLib_CacheParams.hxx>
+
 //! \brief A cache class for Bezier and B-spline curves.
 //!
 //! Defines all data, that can be cached on a span of a curve.
 class BSplCLib_Cache : public Standard_Transient
 {
 public:
-  //! Default constructor
-  Standard_EXPORT BSplCLib_Cache();
-  //! Constructor for caching of 2D curves
+
+  //! Constructor, prepares data structures for caching values on a 2d curve.
   //! \param theDegree     degree of the curve
-  //! \param thePeriodic   identify the curve is periodic
+  //! \param thePeriodic   identify whether the curve is periodic
   //! \param theFlatKnots  knots of Bezier/B-spline curve (with repetitions)
   //! \param thePoles2d    array of poles of 2D curve
   //! \param theWeights    array of weights of corresponding poles
@@ -51,9 +51,10 @@ public:
                                  const TColStd_Array1OfReal&    theFlatKnots,
                                  const TColgp_Array1OfPnt2d&    thePoles2d,
                                  const TColStd_Array1OfReal*    theWeights = NULL);
-  //! Constructor for caching of 3D curves
+
+  //! Constructor, prepares data structures for caching values on a 3d curve.
   //! \param theDegree     degree of the curve
-  //! \param thePeriodic   identify the curve is periodic
+  //! \param thePeriodic   identify whether the curve is periodic
   //! \param theFlatKnots  knots of Bezier/B-spline curve (with repetitions)
   //! \param thePoles      array of poles of 3D curve
   //! \param theWeights    array of weights of corresponding poles
@@ -69,27 +70,20 @@ public:
 
   //! Recomputes the cache data for 2D curves. Does not verify validity of the cache
   //! \param theParameter  the value on the knot's axis to identify the span
-  //! \param theDegree     degree of the curve
-  //! \param thePeriodic   identify the curve is periodic
   //! \param theFlatKnots  knots of Bezier/B-spline curve (with repetitions)
   //! \param thePoles2d    array of poles of 2D curve
   //! \param theWeights    array of weights of corresponding poles
   Standard_EXPORT void BuildCache(const Standard_Real&           theParameter,
-                                  const Standard_Integer&        theDegree,
-                                  const Standard_Boolean&        thePeriodic,
                                   const TColStd_Array1OfReal&    theFlatKnots,
                                   const TColgp_Array1OfPnt2d&    thePoles2d,
-                                  const TColStd_Array1OfReal*    theWeights = NULL);
+                                  const TColStd_Array1OfReal*    theWeights);
+
   //! Recomputes the cache data for 3D curves. Does not verify validity of the cache
   //! \param theParameter  the value on the knot's axis to identify the span
-  //! \param theDegree     degree of the curve
-  //! \param thePeriodic   identify the curve is periodic
   //! \param theFlatKnots  knots of Bezier/B-spline curve (with repetitions)
   //! \param thePoles      array of poles of 3D curve
   //! \param theWeights    array of weights of corresponding poles
   Standard_EXPORT void BuildCache(const Standard_Real&           theParameter,
-                                  const Standard_Integer&        theDegree,
-                                  const Standard_Boolean&        thePeriodic,
                                   const TColStd_Array1OfReal&    theFlatKnots,
                                   const TColgp_Array1OfPnt&      thePoles,
                                   const TColStd_Array1OfReal*    theWeights = NULL);
@@ -142,10 +136,6 @@ public:
   DEFINE_STANDARD_RTTIEXT(BSplCLib_Cache,Standard_Transient)
 
 protected:
-  //! Normalizes the parameter for periodical curves
-  //! \param theFlatKnots knots with repetitions
-  //! \param theParameter the value to be normalized into the knots array
-  void PeriodicNormalization(const TColStd_Array1OfReal& theFlatKnots, Standard_Real& theParameter) const;
 
   //! Fills array of derivatives in the selected point of the curve
   //! \param[in]  theParameter  parameter of the calculation
@@ -156,21 +146,18 @@ protected:
                            const Standard_Integer& theDerivative, 
                                  Standard_Real&    theDerivArray) const;
 
+  // copying is prohibited
+  BSplCLib_Cache (const BSplCLib_Cache&);
+  void operator = (const BSplCLib_Cache&);
+
 private:
-  Handle(TColStd_HArray2OfReal) myPolesWeights; ///< array of poles and weights of calculated cache
+  Standard_Boolean myIsRational;                //!< identifies the rationality of Bezier/B-spline curve
+  BSplCLib_CacheParams myParams;                //!< cache parameters
+  Handle(TColStd_HArray2OfReal) myPolesWeights; //!< array of poles and weights of calculated cache
                                                 // the array has following structure:
                                                 //       x1 y1 [z1] [w1]
                                                 //       x2 y2 [z2] [w2] etc
                                                 // for 2D-curves there is no z conponent, for non-rational curves there is no weight
-
-  Standard_Boolean              myIsRational; ///< identifies the rationality of Bezier/B-spline curve
-  Standard_Real                 mySpanStart;  ///< parameter for the first point of the span
-  Standard_Real                 mySpanLength; ///< length of the span
-  Standard_Integer              mySpanIndex;  ///< index of the span on Bezier/B-spline curve
-  Standard_Integer              mySpanIndexMin; ///< minimal index of span on Bezier/B-spline curve
-  Standard_Integer              mySpanIndexMax; ///< maximal number of spans on Bezier/B-spline curve
-  Standard_Integer              myDegree;     ///< degree of Bezier/B-spline
-  Handle(TColStd_HArray1OfReal) myFlatKnots;  ///< knots of Bezier/B-spline (used for periodic normalization of parameters, exists only for periodical splines)
 };
 
 DEFINE_STANDARD_HANDLE(BSplCLib_Cache, Standard_Transient)
diff --git a/src/BSplCLib/BSplCLib_CacheParams.hxx b/src/BSplCLib/BSplCLib_CacheParams.hxx
new file mode 100644 (file)
index 0000000..a25094b
--- /dev/null
@@ -0,0 +1,106 @@
+// Copyright (c) 2018 OPEN CASCADE SAS
+//
+// This file is part of Open CASCADE Technology software library.
+//
+// This library is free software; you can redistribute it and/or modify it under
+// the terms of the GNU Lesser General Public License version 2.1 as published
+// by the Free Software Foundation, with special exception defined in the file
+// OCCT_LGPL_EXCEPTION.txt. Consult the file LICENSE_LGPL_21.txt included in OCCT
+// distribution for complete text of the license and disclaimer of any warranty.
+//
+// Alternatively, this file may be used under the terms of Open CASCADE
+// commercial license or contractual agreement.
+
+#ifndef _BSplCLib_CacheParams_Headerfile
+#define _BSplCLib_CacheParams_Headerfile
+
+#include <Standard_Real.hxx>
+#include <TColStd_Array1OfReal.hxx>
+
+#include <BSplCLib.hxx>
+
+//! Simple structure containing parameters describing parameterization
+//! of a B-spline curve or a surface in one direction (U or V),
+//! and data of the current span for its caching
+struct BSplCLib_CacheParams
+{
+  const Standard_Integer Degree;      ///< degree of Bezier/B-spline
+  const Standard_Boolean IsPeriodic;  ///< true of the B-spline is periodic
+  const Standard_Real FirstParameter; ///< first valid parameter
+  const Standard_Real LastParameter;  ///< last valid parameter
+
+  const Standard_Integer SpanIndexMin; ///< minimal index of span
+  const Standard_Integer SpanIndexMax; ///< maximal index of span
+
+  Standard_Real    SpanStart;    ///< parameter for the frst point of the span
+  Standard_Real    SpanLength;   ///< length of the span
+  Standard_Integer SpanIndex;    ///< index of the span
+
+  //! Constructor, prepares data structures for caching.
+  //! \param theDegree     degree of the B-spline (or Bezier)
+  //! \param thePeriodic   identify whether the B-spline is periodic
+  //! \param theFlatKnots  knots of Bezier / B-spline parameterization
+  BSplCLib_CacheParams (Standard_Integer theDegree, Standard_Boolean thePeriodic,
+                        const TColStd_Array1OfReal& theFlatKnots)
+  : Degree(theDegree),
+    IsPeriodic(thePeriodic),
+    FirstParameter(theFlatKnots.Value(theFlatKnots.Lower() + theDegree)),
+    LastParameter(theFlatKnots.Value(theFlatKnots.Upper() - theDegree)),
+    SpanIndexMin(theFlatKnots.Lower() + theDegree),
+    SpanIndexMax(theFlatKnots.Upper() - theDegree - 1),
+    SpanStart(0.),
+    SpanLength(0.),
+    SpanIndex(0)
+  {}
+
+  //! Normalizes the parameter for periodic B-splines
+  //! \param theParameter the value to be normalized into the knots array
+  Standard_Real PeriodicNormalization (Standard_Real theParameter) const
+  {
+    if (IsPeriodic)
+    {
+      if (theParameter < FirstParameter)
+      {
+        Standard_Real aPeriod = LastParameter - FirstParameter;
+        Standard_Real aScale = IntegerPart ((FirstParameter - theParameter) / aPeriod);
+        return theParameter + aPeriod * (aScale + 1.0);
+      }
+      if (theParameter > LastParameter)
+      {
+        Standard_Real aPeriod = LastParameter - FirstParameter;
+        Standard_Real aScale = IntegerPart ((theParameter - LastParameter) / aPeriod);
+        return theParameter - aPeriod * (aScale + 1.0);
+      }
+    }
+    return theParameter;
+  }
+
+  //! Verifies validity of the cache using flat parameter of the point
+  //! \param theParameter parameter of the point placed in the span
+  Standard_Boolean IsCacheValid (Standard_Real theParameter) const
+  {
+    Standard_Real aNewParam = PeriodicNormalization  (theParameter);
+    Standard_Real aDelta = aNewParam - SpanStart;
+    return ((aDelta >= 0.0 || SpanIndex == SpanIndexMin) &&
+            (aDelta < SpanLength || SpanIndex == SpanIndexMax));
+  }
+
+  //! Computes span for the specified parameter
+  //! \param theParameter parameter of the point placed in the span
+  //! \param theFlatKnots  knots of Bezier / B-spline parameterization
+  void LocateParameter (Standard_Real& theParameter, const TColStd_Array1OfReal& theFlatKnots)
+  {
+    SpanIndex = 0;
+    BSplCLib::LocateParameter (Degree, theFlatKnots, BSplCLib::NoMults(), 
+                               theParameter, IsPeriodic, SpanIndex, theParameter);
+    SpanStart  = theFlatKnots.Value(SpanIndex);
+    SpanLength = theFlatKnots.Value(SpanIndex + 1) - SpanStart;
+  }
+
+private:
+  // copying is prohibited
+  BSplCLib_CacheParams (const BSplCLib_CacheParams&);
+  void operator = (const BSplCLib_CacheParams&);
+};
+
+#endif
index 2975ef6..5aec2f5 100755 (executable)
@@ -7,6 +7,7 @@ BSplCLib_3.cxx
 BSplCLib_BzSyntaxes.cxx
 BSplCLib_Cache.cxx
 BSplCLib_Cache.hxx
+BSplCLib_CacheParams.hxx
 BSplCLib_CurveComputation.gxx
 BSplCLib_EvaluatorFunction.hxx
 BSplCLib_KnotDistribution.hxx
index dd9ba2e..7ccce9a 100644 (file)
@@ -31,154 +31,60 @@ static Standard_Real* ConvertArray(const Handle(TColStd_HArray2OfReal)& theHArra
   return (Standard_Real*) &(anArray(anArray.LowerRow(), anArray.LowerCol()));
 }
 
-
-BSplSLib_Cache::BSplSLib_Cache()
-{
-  myPolesWeights.Nullify();
-  myIsRational = Standard_False;
-  mySpanStart[0]  = mySpanStart[1]  = 0.0;
-  mySpanLength[0] = mySpanLength[1] = 0.0;
-  mySpanIndex[0]  = mySpanIndex[1]  = 0;
-  myDegree[0]     = myDegree[1]     = 0;
-  myFlatKnots[0].Nullify();
-  myFlatKnots[1].Nullify();
-}
-
 BSplSLib_Cache::BSplSLib_Cache(const Standard_Integer&        theDegreeU,
                                const Standard_Boolean&        thePeriodicU,
                                const TColStd_Array1OfReal&    theFlatKnotsU,
                                const Standard_Integer&        theDegreeV,
                                const Standard_Boolean&        thePeriodicV,
                                const TColStd_Array1OfReal&    theFlatKnotsV,
-                               const TColgp_Array2OfPnt&      thePoles,
                                const TColStd_Array2OfReal*    theWeights)
+: myIsRational(theWeights != NULL),
+  myParamsU (theDegreeU, thePeriodicU, theFlatKnotsU),
+  myParamsV (theDegreeV, thePeriodicV, theFlatKnotsV)
 {
-  Standard_Real aU = theFlatKnotsU.Value(theFlatKnotsU.Lower() + theDegreeU);
-  Standard_Real aV = theFlatKnotsV.Value(theFlatKnotsV.Lower() + theDegreeV);
-
-  BuildCache(aU, aV, 
-             theDegreeU, thePeriodicU, theFlatKnotsU, 
-             theDegreeV, thePeriodicV, theFlatKnotsV, 
-             thePoles, theWeights);
+  Standard_Integer aMinDegree = Min (theDegreeU, theDegreeV);
+  Standard_Integer aMaxDegree = Max (theDegreeU, theDegreeV);
+  Standard_Integer aPWColNumber = (myIsRational ? 4 : 3);
+  myPolesWeights = new TColStd_HArray2OfReal(1, aMaxDegree + 1, 1, aPWColNumber * (aMinDegree + 1));
 }
 
-
 Standard_Boolean BSplSLib_Cache::IsCacheValid(Standard_Real theParameterU,
                                               Standard_Real theParameterV) const
 {
-  Standard_Real aNewU = theParameterU;
-  Standard_Real aNewV = theParameterV;
-  if (!myFlatKnots[0].IsNull())
-    PeriodicNormalization(myDegree[0], myFlatKnots[0]->Array1(), aNewU);
-  if (!myFlatKnots[1].IsNull())
-    PeriodicNormalization(myDegree[1], myFlatKnots[1]->Array1(), aNewV);
-
-  Standard_Real aDelta0 = aNewU - mySpanStart[0];
-  Standard_Real aDelta1 = aNewV - mySpanStart[1];
-  return ((aDelta0 >= -mySpanLength[0] || mySpanIndex[0] == mySpanIndexMin[0]) &&
-          (aDelta0 < mySpanLength[0] || mySpanIndex[0] == mySpanIndexMax[0]) &&
-          (aDelta1 >= -mySpanLength[1] || mySpanIndex[1] == mySpanIndexMin[1]) &&
-          (aDelta1 < mySpanLength[1] || mySpanIndex[1] == mySpanIndexMax[1]));
-}
-
-void BSplSLib_Cache::PeriodicNormalization(const Standard_Integer& theDegree, 
-                                           const TColStd_Array1OfReal& theFlatKnots, 
-                                           Standard_Real& theParameter) const
-{
-  Standard_Real aPeriod = theFlatKnots.Value(theFlatKnots.Upper() - theDegree) - 
-                          theFlatKnots.Value(theDegree + 1) ;
-  if (theParameter < theFlatKnots.Value(theDegree + 1))
-  {
-    Standard_Real aScale = IntegerPart(
-        (theFlatKnots.Value(theDegree + 1) - theParameter) / aPeriod);
-    theParameter += aPeriod * (aScale + 1.0);
-  }
-  if (theParameter > theFlatKnots.Value(theFlatKnots.Upper() - theDegree))
-  {
-    Standard_Real aScale = IntegerPart(
-        (theParameter - theFlatKnots.Value(theFlatKnots.Upper() - theDegree)) / aPeriod);
-    theParameter -= aPeriod * (aScale + 1.0);
-  }
+  return myParamsU.IsCacheValid (theParameterU) && 
+         myParamsV.IsCacheValid (theParameterV);
 }
 
-
-void BSplSLib_Cache::BuildCache(const Standard_Real&           theParameterU, 
-                                const Standard_Real&           theParameterV, 
-                                const Standard_Integer&        theDegreeU, 
-                                const Standard_Boolean&        thePeriodicU, 
-                                const TColStd_Array1OfReal&    theFlatKnotsU, 
-                                const Standard_Integer&        theDegreeV, 
-                                const Standard_Boolean&        thePeriodicV, 
-                                const TColStd_Array1OfReal&    theFlatKnotsV, 
-                                const TColgp_Array2OfPnt&      thePoles, 
+void BSplSLib_Cache::BuildCache(const Standard_Real&           theParameterU,
+                                const Standard_Real&           theParameterV,
+                                const TColStd_Array1OfReal&    theFlatKnotsU,
+                                const TColStd_Array1OfReal&    theFlatKnotsV,
+                                const TColgp_Array2OfPnt&      thePoles,
                                 const TColStd_Array2OfReal*    theWeights)
 {
   // Normalize the parameters for periodical B-splines
-  Standard_Real aNewParamU = theParameterU;
-  if (thePeriodicU)
-  {
-    PeriodicNormalization(theDegreeU, theFlatKnotsU, aNewParamU);
-    myFlatKnots[0] = new TColStd_HArray1OfReal(1, theFlatKnotsU.Length());
-    myFlatKnots[0]->ChangeArray1() = theFlatKnotsU;
-  }
-  else if (!myFlatKnots[0].IsNull()) // Periodical curve became non-periodical
-    myFlatKnots[0].Nullify();
-
-  Standard_Real aNewParamV = theParameterV;
-  if (thePeriodicV)
-  {
-    PeriodicNormalization(theDegreeV, theFlatKnotsV, aNewParamV);
-    myFlatKnots[1] = new TColStd_HArray1OfReal(1, theFlatKnotsV.Length());
-    myFlatKnots[1]->ChangeArray1() = theFlatKnotsV;
-  }
-  else if (!myFlatKnots[1].IsNull()) // Periodical curve became non-periodical
-    myFlatKnots[1].Nullify();
-
-  Standard_Integer aMinDegree = Min(theDegreeU, theDegreeV);
-  Standard_Integer aMaxDegree = Max(theDegreeU, theDegreeV);
-
-  // Change the size of cached data if needed
-  myIsRational = (theWeights != NULL);
-  Standard_Integer aPWColNumber = myIsRational ? 4 : 3;
-  if (theDegreeU > myDegree[0] || theDegreeV > myDegree[1])
-    myPolesWeights = new TColStd_HArray2OfReal(1, aMaxDegree + 1, 1, aPWColNumber * (aMinDegree + 1));
-
-  myDegree[0] = theDegreeU;
-  myDegree[1] = theDegreeV;
-  mySpanIndex[0] = mySpanIndex[1] = 0;
-  BSplCLib::LocateParameter(theDegreeU, theFlatKnotsU, BSplCLib::NoMults(), aNewParamU, 
-                            thePeriodicU, mySpanIndex[0], aNewParamU);
-  BSplCLib::LocateParameter(theDegreeV, theFlatKnotsV, BSplCLib::NoMults(), aNewParamV, 
-                            thePeriodicV, mySpanIndex[1], aNewParamV);
-
-  // Protection against Out of Range exception.
-  if (mySpanIndex[0] >= theFlatKnotsU.Length()) {
-    mySpanIndex[0] = theFlatKnotsU.Length() - 1;
-  }
-
-  mySpanLength[0] = (theFlatKnotsU.Value(mySpanIndex[0] + 1) - theFlatKnotsU.Value(mySpanIndex[0])) * 0.5;
-  mySpanStart[0]  = theFlatKnotsU.Value(mySpanIndex[0]) + mySpanLength[0];
+  Standard_Real aNewParamU = myParamsU.PeriodicNormalization (theParameterU);
+  Standard_Real aNewParamV = myParamsV.PeriodicNormalization (theParameterV);
 
-  // Protection against Out of Range exception.
-  if (mySpanIndex[1] >= theFlatKnotsV.Length()) {
-    mySpanIndex[1] = theFlatKnotsV.Length() - 1;
-  }
+  myParamsU.LocateParameter (aNewParamU, theFlatKnotsU);
+  myParamsV.LocateParameter (aNewParamV, theFlatKnotsV);
 
-  mySpanLength[1] = (theFlatKnotsV.Value(mySpanIndex[1] + 1) - theFlatKnotsV.Value(mySpanIndex[1])) * 0.5;
-  mySpanStart[1]  = theFlatKnotsV.Value(mySpanIndex[1]) + mySpanLength[1];
-  mySpanIndexMin[0] = thePeriodicU ? 0 : theDegreeU + 1;
-  mySpanIndexMax[0] = theFlatKnotsU.Length() - 1 - theDegreeU;
-  mySpanIndexMin[1] = thePeriodicV ? 0 : theDegreeV + 1;
-  mySpanIndexMax[1] = theFlatKnotsV.Length() - 1 - theDegreeV;
+  // BSplSLib uses different convention for span parameters than BSplCLib
+  // (Start is in the middle of the span and length is half-span),
+  // thus we need to amend them here
+  Standard_Real aSpanLengthU = 0.5 * myParamsU.SpanLength;
+  Standard_Real aSpanStartU = myParamsU.SpanStart + aSpanLengthU;
+  Standard_Real aSpanLengthV = 0.5 * myParamsV.SpanLength;
+  Standard_Real aSpanStartV = myParamsV.SpanStart + aSpanLengthV;
 
   // Calculate new cache data
-  BSplSLib::BuildCache(mySpanStart[0],  mySpanStart[1], 
-                       mySpanLength[0], mySpanLength[1], 
-                       thePeriodicU,    thePeriodicV, 
-                       theDegreeU,      theDegreeV, 
-                       mySpanIndex[0],  mySpanIndex[1], 
-                       theFlatKnotsU,   theFlatKnotsV, 
-                       thePoles, theWeights, myPolesWeights->ChangeArray2());
+  BSplSLib::BuildCache (aSpanStartU,  aSpanStartV,
+                        aSpanLengthU, aSpanLengthV,
+                        myParamsU.IsPeriodic, myParamsV.IsPeriodic,
+                        myParamsU.Degree,     myParamsV.Degree,
+                        myParamsU.SpanIndex,  myParamsV.SpanIndex,
+                        theFlatKnotsU,        theFlatKnotsV,
+                        thePoles, theWeights, myPolesWeights->ChangeArray2());
 }
 
 
@@ -186,24 +92,29 @@ void BSplSLib_Cache::D0(const Standard_Real& theU,
                         const Standard_Real& theV, 
                               gp_Pnt&        thePoint) const
 {
-  Standard_Real aNewU = theU;
-  Standard_Real aNewV = theV;
-  if (!myFlatKnots[0].IsNull()) // B-spline is U-periodical
-    PeriodicNormalization(myDegree[0], myFlatKnots[0]->Array1(), aNewU);
-  aNewU = (aNewU - mySpanStart[0]) / mySpanLength[0];
-  if (!myFlatKnots[1].IsNull()) // B-spline is V-periodical
-    PeriodicNormalization(myDegree[1], myFlatKnots[1]->Array1(), aNewV);
-  aNewV = (aNewV - mySpanStart[1]) / mySpanLength[1];
+  Standard_Real aNewU = myParamsU.PeriodicNormalization (theU);
+  Standard_Real aNewV = myParamsV.PeriodicNormalization (theV);
+
+  // BSplSLib uses different convention for span parameters than BSplCLib
+  // (Start is in the middle of the span and length is half-span),
+  // thus we need to amend them here
+  Standard_Real aSpanLengthU = 0.5 * myParamsU.SpanLength;
+  Standard_Real aSpanStartU = myParamsU.SpanStart + aSpanLengthU;
+  Standard_Real aSpanLengthV = 0.5 * myParamsV.SpanLength;
+  Standard_Real aSpanStartV = myParamsV.SpanStart + aSpanLengthV;
+
+  aNewU = (aNewU - aSpanStartU) / aSpanLengthU;
+  aNewV = (aNewV - aSpanStartV) / aSpanLengthV;
 
   Standard_Real* aPolesArray = ConvertArray(myPolesWeights);
   Standard_Real aPoint[4];
 
   Standard_Integer aDimension = myIsRational ? 4 : 3;
   Standard_Integer aCacheCols = myPolesWeights->RowLength();
-  Standard_Integer aMinMaxDegree[2] = {Min(myDegree[0], myDegree[1]), 
-                                       Max(myDegree[0], myDegree[1])};
+  Standard_Integer aMinMaxDegree[2] = {Min(myParamsU.Degree, myParamsV.Degree),
+                                       Max(myParamsU.Degree, myParamsV.Degree)};
   Standard_Real aParameters[2];
-  if (myDegree[0] > myDegree[1])
+  if (myParamsU.Degree > myParamsV.Degree)
   {
     aParameters[0] = aNewV;
     aParameters[1] = aNewU;
@@ -238,16 +149,21 @@ void BSplSLib_Cache::D1(const Standard_Real& theU,
                               gp_Vec&        theTangentU, 
                               gp_Vec&        theTangentV) const
 {
-  Standard_Real aNewU = theU;
-  Standard_Real aNewV = theV;
-  Standard_Real anInvU = 1.0 / mySpanLength[0];
-  Standard_Real anInvV = 1.0 / mySpanLength[1];
-  if (!myFlatKnots[0].IsNull()) // B-spline is U-periodical
-    PeriodicNormalization(myDegree[0], myFlatKnots[0]->Array1(), aNewU);
-  aNewU = (aNewU - mySpanStart[0]) * anInvU;
-  if (!myFlatKnots[1].IsNull()) // B-spline is V-periodical
-    PeriodicNormalization(myDegree[1], myFlatKnots[1]->Array1(), aNewV);
-  aNewV = (aNewV - mySpanStart[1]) * anInvV;
+  Standard_Real aNewU = myParamsU.PeriodicNormalization (theU);
+  Standard_Real aNewV = myParamsV.PeriodicNormalization (theV);
+
+  // BSplSLib uses different convention for span parameters than BSplCLib
+  // (Start is in the middle of the span and length is half-span),
+  // thus we need to amend them here
+  Standard_Real aSpanLengthU = 0.5 * myParamsU.SpanLength;
+  Standard_Real aSpanStartU = myParamsU.SpanStart + aSpanLengthU;
+  Standard_Real aSpanLengthV = 0.5 * myParamsV.SpanLength;
+  Standard_Real aSpanStartV = myParamsV.SpanStart + aSpanLengthV;
+
+  Standard_Real anInvU = 1.0 / aSpanLengthU;
+  Standard_Real anInvV = 1.0 / aSpanLengthV;
+  aNewU = (aNewU - aSpanStartU) * anInvU;
+  aNewV = (aNewV - aSpanStartV) * anInvV;
 
   Standard_Real* aPolesArray = ConvertArray(myPolesWeights);
   Standard_Real aPntDeriv[16]; // result storage (point and derivative coordinates)
@@ -255,11 +171,11 @@ void BSplSLib_Cache::D1(const Standard_Real& theU,
 
   Standard_Integer aDimension = myIsRational ? 4 : 3;
   Standard_Integer aCacheCols = myPolesWeights->RowLength();
-  Standard_Integer aMinMaxDegree[2] = {Min(myDegree[0], myDegree[1]), 
-                                       Max(myDegree[0], myDegree[1])};
+  Standard_Integer aMinMaxDegree[2] = {Min(myParamsU.Degree, myParamsV.Degree),
+                                       Max(myParamsU.Degree, myParamsV.Degree)};
 
   Standard_Real aParameters[2];
-  if (myDegree[0] > myDegree[1])
+  if (myParamsU.Degree > myParamsV.Degree)
   {
     aParameters[0] = aNewV;
     aParameters[1] = aNewU;
@@ -293,7 +209,7 @@ void BSplSLib_Cache::D1(const Standard_Real& theU,
   }
 
   thePoint.SetCoord(aResult[0], aResult[1], aResult[2]);
-  if (myDegree[0] > myDegree[1])
+  if (myParamsU.Degree > myParamsV.Degree)
   {
     theTangentV.SetCoord(aResult[aDimension], aResult[aDimension + 1], aResult[aDimension + 2]);
     Standard_Integer aShift = aDimension<<1;
@@ -319,16 +235,21 @@ void BSplSLib_Cache::D2(const Standard_Real& theU,
                               gp_Vec&        theCurvatureV, 
                               gp_Vec&        theCurvatureUV) const
 {
-  Standard_Real aNewU = theU;
-  Standard_Real aNewV = theV;
-  Standard_Real anInvU = 1.0 / mySpanLength[0];
-  Standard_Real anInvV = 1.0 / mySpanLength[1];
-  if (!myFlatKnots[0].IsNull()) // B-spline is U-periodical
-    PeriodicNormalization(myDegree[0], myFlatKnots[0]->Array1(), aNewU);
-  aNewU = (aNewU - mySpanStart[0]) * anInvU;
-  if (!myFlatKnots[1].IsNull()) // B-spline is V-periodical
-    PeriodicNormalization(myDegree[1], myFlatKnots[1]->Array1(), aNewV);
-  aNewV = (aNewV - mySpanStart[1]) * anInvV;
+  Standard_Real aNewU = myParamsU.PeriodicNormalization (theU);
+  Standard_Real aNewV = myParamsV.PeriodicNormalization (theV);
+
+  // BSplSLib uses different convention for span parameters than BSplCLib
+  // (Start is in the middle of the span and length is half-span),
+  // thus we need to amend them here
+  Standard_Real aSpanLengthU = 0.5 * myParamsU.SpanLength;
+  Standard_Real aSpanStartU = myParamsU.SpanStart + aSpanLengthU;
+  Standard_Real aSpanLengthV = 0.5 * myParamsV.SpanLength;
+  Standard_Real aSpanStartV = myParamsV.SpanStart + aSpanLengthV;
+
+  Standard_Real anInvU = 1.0 / aSpanLengthU;
+  Standard_Real anInvV = 1.0 / aSpanLengthV;
+  aNewU = (aNewU - aSpanStartU) * anInvU;
+  aNewV = (aNewV - aSpanStartV) * anInvV;
 
   Standard_Real* aPolesArray = ConvertArray(myPolesWeights);
   Standard_Real aPntDeriv[36]; // result storage (point and derivative coordinates)
@@ -336,11 +257,11 @@ void BSplSLib_Cache::D2(const Standard_Real& theU,
 
   Standard_Integer aDimension = myIsRational ? 4 : 3;
   Standard_Integer aCacheCols = myPolesWeights->RowLength();
-  Standard_Integer aMinMaxDegree[2] = {Min(myDegree[0], myDegree[1]), 
-                                       Max(myDegree[0], myDegree[1])};
+  Standard_Integer aMinMaxDegree[2] = {Min(myParamsU.Degree, myParamsV.Degree),
+                                       Max(myParamsU.Degree, myParamsV.Degree)};
 
   Standard_Real aParameters[2];
-  if (myDegree[0] > myDegree[1])
+  if (myParamsU.Degree > myParamsV.Degree)
   {
     aParameters[0] = aNewV;
     aParameters[1] = aNewU;
@@ -390,7 +311,7 @@ void BSplSLib_Cache::D2(const Standard_Real& theU,
   }
 
   thePoint.SetCoord(aResult[0], aResult[1], aResult[2]);
-  if (myDegree[0] > myDegree[1])
+  if (myParamsU.Degree > myParamsV.Degree)
   {
     theTangentV.SetCoord(aResult[aDimension], aResult[aDimension + 1], aResult[aDimension + 2]);
     Standard_Integer aShift = aDimension<<1;
index 821cf7e..965ba6e 100644 (file)
@@ -29,6 +29,8 @@
 #include <TColStd_Array1OfReal.hxx>
 #include <TColStd_Array2OfReal.hxx>
 
+#include <BSplCLib_CacheParams.hxx>
+
 //! \brief A cache class for Bezier and B-spline surfaces.
 //!
 //! Defines all data, that can be cached on a span of the surface.
@@ -36,8 +38,7 @@
 class BSplSLib_Cache : public Standard_Transient
 {
 public:
-  //! Default constructor
-  Standard_EXPORT BSplSLib_Cache();
+
   //! Constructor for caching of the span for the surface
   //! \param theDegreeU    degree along the first parameter (U) of the surface
   //! \param thePeriodicU  identify the surface is periodical along U axis
@@ -45,7 +46,6 @@ public:
   //! \param theDegreeV    degree alogn the second parameter (V) of the surface
   //! \param thePeriodicV  identify the surface is periodical along V axis
   //! \param theFlatKnotsV knots of the surface (with repetition) along V axis
-  //! \param thePoles      array of poles of the surface
   //! \param theWeights    array of weights of corresponding poles
   Standard_EXPORT BSplSLib_Cache(const Standard_Integer&        theDegreeU,
                                  const Standard_Boolean&        thePeriodicU,
@@ -53,7 +53,6 @@ public:
                                  const Standard_Integer&        theDegreeV,
                                  const Standard_Boolean&        thePeriodicV,
                                  const TColStd_Array1OfReal&    theFlatKnotsV,
-                                 const TColgp_Array2OfPnt&      thePoles,
                                  const TColStd_Array2OfReal*    theWeights = NULL);
 
   //! Verifies validity of the cache using parameters of the point
@@ -73,14 +72,10 @@ public:
   //! \param theFlatKnotsV  flat knots of the surface along V axis
   //! \param thePoles       array of poles of the surface
   //! \param theWeights     array of weights of corresponding poles
-  Standard_EXPORT void BuildCache(const Standard_Real&           theParameterU, 
-                                  const Standard_Real&           theParameterV, 
-                                  const Standard_Integer&        theDegreeU, 
-                                  const Standard_Boolean&        thePeriodicU, 
-                                  const TColStd_Array1OfReal&    theFlatKnotsU, 
-                                  const Standard_Integer&        theDegreeV, 
-                                  const Standard_Boolean&        thePeriodicV, 
-                                  const TColStd_Array1OfReal&    theFlatKnotsV, 
+  Standard_EXPORT void BuildCache(const Standard_Real&           theParameterU,
+                                  const Standard_Real&           theParameterV,
+                                  const TColStd_Array1OfReal&    theFlatKnotsU,
+                                  const TColStd_Array1OfReal&    theFlatKnotsV,
                                   const TColgp_Array2OfPnt&      thePoles, 
                                   const TColStd_Array2OfReal*    theWeights = NULL);
 
@@ -123,32 +118,20 @@ public:
 
   DEFINE_STANDARD_RTTIEXT(BSplSLib_Cache,Standard_Transient)
 
-protected:
-  //! Normalizes the parameter for periodical surfaces
-  //! \param[in]     theDegree     degree along selected direction
-  //! \param[in]     theFlatKnots  knots with repetitions along selected direction
-  //! \param[in,out] theParameter  the value to be normalized into the knots array
-  void PeriodicNormalization(const Standard_Integer& theDegree, 
-                             const TColStd_Array1OfReal& theFlatKnots, 
-                                   Standard_Real& theParameter) const;
+private:
+  // copying is prohibited
+  BSplSLib_Cache (const BSplSLib_Cache&);
+  void operator = (const BSplSLib_Cache&);
 
 private:
-  Handle(TColStd_HArray2OfReal) myPolesWeights; ///< array of poles and weights of calculated cache
+  Standard_Boolean myIsRational;                //!< identifies the rationality of Bezier/B-spline surface
+  BSplCLib_CacheParams myParamsU, myParamsV;    //!< cach parameters by U and V directions
+  Handle(TColStd_HArray2OfReal) myPolesWeights; //!< array of poles and weights of calculated cache
                                                 // the array has following structure:
                                                 //       x11 y11 z11 [w11] x12 y12 z12 [w12] ...
                                                 //       x21 y21 z21 [w21] x22 y22 z22 [w22] etc
                                                 // for non-rational surfaces there is no weight;
                                                 // size of array: (max(myDegree)+1) * A*(min(myDegree)+1), where A = 4 or 3
-
-  Standard_Boolean              myIsRational;    ///< identifies the rationality of Bezier/B-spline surface
-  Standard_Real                 mySpanStart[2];  ///< parameters (u, v) for the frst point of the span
-  Standard_Real                 mySpanLength[2]; ///< lengths of the span along corresponding parameter
-  Standard_Integer              mySpanIndex[2];  ///< indexes of the span on Bezier/B-spline surface
-  Standard_Integer              mySpanIndexMin[2]; ///< minimal indexes of span
-  Standard_Integer              mySpanIndexMax[2]; ///< maximal indexes of span
-  Standard_Integer              myDegree[2];     ///< degrees of Bezier/B-spline for each parameter
-  Handle(TColStd_HArray1OfReal) myFlatKnots[2];  ///< arrays of knots of Bezier/B-spline 
-                                                 // (used for periodic normalization of parameters, Null for non-periodical splines)
 };
 
 DEFINE_STANDARD_HANDLE(BSplSLib_Cache, Standard_Transient)
index e4a3074..ad98572 100644 (file)
@@ -556,20 +556,18 @@ void Geom2dAdaptor_Curve::RebuildCache(const Standard_Real theParameter) const
     Standard_Integer aDeg = aBezier->Degree();
     TColStd_Array1OfReal aFlatKnots(BSplCLib::FlatBezierKnots(aDeg), 1, 2 * (aDeg + 1));
     if (myCurveCache.IsNull())
-      myCurveCache = new BSplCLib_Cache(aDeg, aBezier->IsPeriodic(), aFlatKnots,
-        aBezier->Poles(), aBezier->Weights());
-    myCurveCache->BuildCache(theParameter, aDeg, aBezier->IsPeriodic(), aFlatKnots,
-      aBezier->Poles(), aBezier->Weights());
+      myCurveCache = new BSplCLib_Cache (aDeg, aBezier->IsPeriodic(), aFlatKnots,
+                                         aBezier->Poles(), aBezier->Weights());
+    myCurveCache->BuildCache (theParameter, aFlatKnots, aBezier->Poles(), aBezier->Weights());
   }
   else if (myTypeCurve == GeomAbs_BSplineCurve)
   {
     // Create cache for B-spline
     if (myCurveCache.IsNull())
-      myCurveCache = new BSplCLib_Cache(myBSplineCurve->Degree(), myBSplineCurve->IsPeriodic(),
+      myCurveCache = new BSplCLib_Cache (myBSplineCurve->Degree(), myBSplineCurve->IsPeriodic(),
         myBSplineCurve->KnotSequence(), myBSplineCurve->Poles(), myBSplineCurve->Weights());
-    myCurveCache->BuildCache(theParameter, myBSplineCurve->Degree(),
-      myBSplineCurve->IsPeriodic(), myBSplineCurve->KnotSequence(),
-      myBSplineCurve->Poles(), myBSplineCurve->Weights());
+    myCurveCache->BuildCache (theParameter, myBSplineCurve->KnotSequence(),
+                              myBSplineCurve->Poles(), myBSplineCurve->Weights());
   }
 }
 
index 42bf072..0215eef 100644 (file)
@@ -539,8 +539,7 @@ void GeomAdaptor_Curve::RebuildCache(const Standard_Real theParameter) const
     if (myCurveCache.IsNull())
       myCurveCache = new BSplCLib_Cache(aDeg, aBezier->IsPeriodic(), aFlatKnots,
         aBezier->Poles(), aBezier->Weights());
-    myCurveCache->BuildCache(theParameter, aDeg, aBezier->IsPeriodic(), aFlatKnots,
-        aBezier->Poles(), aBezier->Weights());
+    myCurveCache->BuildCache (theParameter, aFlatKnots, aBezier->Poles(), aBezier->Weights());
 }
   else if (myTypeCurve == GeomAbs_BSplineCurve)
 {
@@ -548,9 +547,8 @@ void GeomAdaptor_Curve::RebuildCache(const Standard_Real theParameter) const
     if (myCurveCache.IsNull())
       myCurveCache = new BSplCLib_Cache(myBSplineCurve->Degree(), myBSplineCurve->IsPeriodic(),
         myBSplineCurve->KnotSequence(), myBSplineCurve->Poles(), myBSplineCurve->Weights());
-    myCurveCache->BuildCache(theParameter, myBSplineCurve->Degree(),
-        myBSplineCurve->IsPeriodic(), myBSplineCurve->KnotSequence(),
-        myBSplineCurve->Poles(), myBSplineCurve->Weights());
+    myCurveCache->BuildCache (theParameter, myBSplineCurve->KnotSequence(),
+                              myBSplineCurve->Poles(), myBSplineCurve->Weights());
 }
 }
 
index 04a7b35..ff860f6 100644 (file)
@@ -670,12 +670,9 @@ void GeomAdaptor_Surface::RebuildCache(const Standard_Real theU,
     if (mySurfaceCache.IsNull())
       mySurfaceCache = new BSplSLib_Cache(
         aDegU, aBezier->IsUPeriodic(), aFlatKnotsU,
-        aDegV, aBezier->IsVPeriodic(), aFlatKnotsV,
-        aBezier->Poles(), aBezier->Weights());
-    mySurfaceCache->BuildCache(theU, theV,
-      aDegU, aBezier->IsUPeriodic(), aFlatKnotsU,
-      aDegV, aBezier->IsVPeriodic(), aFlatKnotsV,
-      aBezier->Poles(), aBezier->Weights());
+        aDegV, aBezier->IsVPeriodic(), aFlatKnotsV, aBezier->Weights());
+    mySurfaceCache->BuildCache (theU, theV, aFlatKnotsU, aFlatKnotsV,
+                                aBezier->Poles(), aBezier->Weights());
   }
   else if (mySurfaceType == GeomAbs_BSplineSurface)
   {
@@ -684,11 +681,9 @@ void GeomAdaptor_Surface::RebuildCache(const Standard_Real theU,
       mySurfaceCache = new BSplSLib_Cache(
         myBSplineSurface->UDegree(), myBSplineSurface->IsUPeriodic(), myBSplineSurface->UKnotSequence(),
         myBSplineSurface->VDegree(), myBSplineSurface->IsVPeriodic(), myBSplineSurface->VKnotSequence(),
-        myBSplineSurface->Poles(), myBSplineSurface->Weights());
-    mySurfaceCache->BuildCache(theU, theV,
-      myBSplineSurface->UDegree(), myBSplineSurface->IsUPeriodic(), myBSplineSurface->UKnotSequence(),
-      myBSplineSurface->VDegree(), myBSplineSurface->IsVPeriodic(), myBSplineSurface->VKnotSequence(),
-      myBSplineSurface->Poles(), myBSplineSurface->Weights());
+        myBSplineSurface->Weights());
+    mySurfaceCache->BuildCache (theU, theV, myBSplineSurface->UKnotSequence(), myBSplineSurface->VKnotSequence(),
+                                myBSplineSurface->Poles(), myBSplineSurface->Weights());
   }
 }