1 // Copyright (c) 2014 OPEN CASCADE SAS
3 // This file is part of Open CASCADE Technology software library.
5 // This library is free software; you can redistribute it and/or modify it under
6 // the terms of the GNU Lesser General Public License version 2.1 as published
7 // by the Free Software Foundation, with special exception defined in the file
8 // OCCT_LGPL_EXCEPTION.txt. Consult the file LICENSE_LGPL_21.txt included in OCCT
9 // distribution for complete text of the license and disclaimer of any warranty.
11 // Alternatively, this file may be used under the terms of Open CASCADE
12 // commercial license or contractual agreement.
14 #include <BSplCLib_Cache.hxx>
15 #include <BSplCLib.hxx>
17 #include <NCollection_LocalArray.hxx>
19 #include <TColgp_HArray1OfPnt.hxx>
20 #include <TColgp_HArray1OfPnt2d.hxx>
21 #include <TColStd_HArray1OfReal.hxx>
22 #include <TColStd_HArray2OfReal.hxx>
25 //! Converts handle of array of Standard_Real into the pointer to Standard_Real
26 static Standard_Real* ConvertArray(const Handle(TColStd_HArray2OfReal)& theHArray)
28 const TColStd_Array2OfReal& anArray = theHArray->Array2();
29 return (Standard_Real*) &(anArray(anArray.LowerRow(), anArray.LowerCol()));
33 BSplCLib_Cache::BSplCLib_Cache()
35 myPolesWeights.Nullify();
36 myIsRational = Standard_False;
41 myFlatKnots.Nullify();
44 BSplCLib_Cache::BSplCLib_Cache(const Standard_Integer& theDegree,
45 const Standard_Boolean& thePeriodic,
46 const TColStd_Array1OfReal& theFlatKnots,
47 const TColgp_Array1OfPnt2d& thePoles2d,
48 const TColStd_Array1OfReal& theWeights)
50 Standard_Real aCacheParam = theFlatKnots.Value(theFlatKnots.Lower() + theDegree);
51 BuildCache(aCacheParam, theDegree, thePeriodic,
52 theFlatKnots, thePoles2d, theWeights);
55 BSplCLib_Cache::BSplCLib_Cache(const Standard_Integer& theDegree,
56 const Standard_Boolean& thePeriodic,
57 const TColStd_Array1OfReal& theFlatKnots,
58 const TColgp_Array1OfPnt& thePoles,
59 const TColStd_Array1OfReal& theWeights)
61 Standard_Real aCacheParam = theFlatKnots.Value(theFlatKnots.Lower() + theDegree);
62 BuildCache(aCacheParam, theDegree, thePeriodic,
63 theFlatKnots, thePoles, theWeights);
67 Standard_Boolean BSplCLib_Cache::IsCacheValid(Standard_Real theParameter) const
69 Standard_Real aNewParam = theParameter;
70 if (!myFlatKnots.IsNull())
71 PeriodicNormalization(myFlatKnots->Array1(), aNewParam);
73 Standard_Real aDelta = aNewParam - mySpanStart;
74 return (aDelta >= 0.0 && (aDelta < mySpanLength || mySpanIndex == mySpanIndexMax));
77 void BSplCLib_Cache::PeriodicNormalization(const TColStd_Array1OfReal& theFlatKnots,
78 Standard_Real& theParameter) const
80 Standard_Real aPeriod = theFlatKnots.Value(theFlatKnots.Upper() - myDegree) -
81 theFlatKnots.Value(myDegree + 1) ;
82 if (theParameter < theFlatKnots.Value(myDegree + 1))
84 Standard_Real aScale = IntegerPart(
85 (theFlatKnots.Value(myDegree + 1) - theParameter) / aPeriod);
86 theParameter += aPeriod * (aScale + 1.0);
88 if (theParameter > theFlatKnots.Value(theFlatKnots.Upper() - myDegree))
90 Standard_Real aScale = IntegerPart(
91 (theParameter - theFlatKnots.Value(theFlatKnots.Upper() - myDegree)) / aPeriod);
92 theParameter -= aPeriod * (aScale + 1.0);
97 void BSplCLib_Cache::BuildCache(const Standard_Real& theParameter,
98 const Standard_Integer& theDegree,
99 const Standard_Boolean& thePeriodic,
100 const TColStd_Array1OfReal& theFlatKnots,
101 const TColgp_Array1OfPnt2d& thePoles2d,
102 const TColStd_Array1OfReal& theWeights)
104 // Normalize theParameter for periodical B-splines
105 Standard_Real aNewParam = theParameter;
108 PeriodicNormalization(theFlatKnots, aNewParam);
109 myFlatKnots = new TColStd_HArray1OfReal(1, theFlatKnots.Length());
110 myFlatKnots->ChangeArray1() = theFlatKnots;
112 else if (!myFlatKnots.IsNull()) // Periodical curve became non-periodical
113 myFlatKnots.Nullify();
115 // Change the size of cached data if needed
116 myIsRational = (&theWeights != NULL);
117 Standard_Integer aPWColNumber = myIsRational ? 3 : 2;
118 if (theDegree > myDegree)
119 myPolesWeights = new TColStd_HArray2OfReal(1, theDegree + 1, 1, aPWColNumber);
121 myDegree = theDegree;
123 BSplCLib::LocateParameter(theDegree, theFlatKnots, BSplCLib::NoMults(),
124 aNewParam, thePeriodic, mySpanIndex, aNewParam);
125 mySpanStart = theFlatKnots.Value(mySpanIndex);
126 mySpanLength = theFlatKnots.Value(mySpanIndex + 1) - mySpanStart;
127 mySpanIndexMax = theFlatKnots.Length() - 1 - theDegree;
129 // Calculate new cache data
130 BSplCLib::BuildCache(mySpanStart, mySpanLength, thePeriodic, theDegree,
131 theFlatKnots, thePoles2d, theWeights,
132 myPolesWeights->ChangeArray2());
135 void BSplCLib_Cache::BuildCache(const Standard_Real& theParameter,
136 const Standard_Integer& theDegree,
137 const Standard_Boolean& thePeriodic,
138 const TColStd_Array1OfReal& theFlatKnots,
139 const TColgp_Array1OfPnt& thePoles,
140 const TColStd_Array1OfReal& theWeights)
142 // Create list of knots with repetitions and normalize theParameter for periodical B-splines
143 Standard_Real aNewParam = theParameter;
146 PeriodicNormalization(theFlatKnots, aNewParam);
147 myFlatKnots = new TColStd_HArray1OfReal(1, theFlatKnots.Length());
148 myFlatKnots->ChangeArray1() = theFlatKnots;
150 else if (!myFlatKnots.IsNull()) // Periodical curve became non-periodical
151 myFlatKnots.Nullify();
153 // Change the size of cached data if needed
154 myIsRational = (&theWeights != NULL);
155 Standard_Integer aPWColNumber = myIsRational ? 4 : 3;
156 if (theDegree > myDegree)
157 myPolesWeights = new TColStd_HArray2OfReal(1, theDegree + 1, 1, aPWColNumber);
159 myDegree = theDegree;
161 BSplCLib::LocateParameter(theDegree, theFlatKnots, BSplCLib::NoMults(),
162 aNewParam, thePeriodic, mySpanIndex, aNewParam);
163 mySpanStart = theFlatKnots.Value(mySpanIndex);
164 mySpanLength = theFlatKnots.Value(mySpanIndex + 1) - mySpanStart;
165 mySpanIndexMax = theFlatKnots.Length() - 1 - theDegree;
167 // Calculate new cache data
168 BSplCLib::BuildCache(mySpanStart, mySpanLength, thePeriodic, theDegree,
169 theFlatKnots, thePoles, theWeights,
170 myPolesWeights->ChangeArray2());
174 void BSplCLib_Cache::CalculateDerivative(const Standard_Real& theParameter,
175 const Standard_Integer& theDerivative,
176 Standard_Real& theDerivArray) const
178 Standard_Real aNewParameter = theParameter;
179 if (!myFlatKnots.IsNull()) // B-spline is periodical
180 PeriodicNormalization(myFlatKnots->Array1(), aNewParameter);
181 aNewParameter = (aNewParameter - mySpanStart) / mySpanLength;
183 Standard_Real* aPolesArray = ConvertArray(myPolesWeights);
184 Standard_Integer aDimension = myPolesWeights->RowLength(); // number of columns
186 // Temporary container. The maximal size of this container is defined by:
187 // 1) maximal derivative for cache evaluation, which is 3, plus one row for function values,
188 // 2) and maximal dimension of the point, which is 3, plus one column for weights.
189 Standard_Real aTmpContainer[16];
191 // When the PLib::RationaDerivative needs to be called, use temporary container
192 Standard_Real* aPntDeriv = myIsRational ? aTmpContainer : &theDerivArray;
194 // When the degree of curve is lesser than the requested derivative,
195 // nullify array cells corresponding to greater derivatives
196 Standard_Integer aDerivative = theDerivative;
197 if (myDegree < theDerivative)
199 aDerivative = myDegree;
200 for (Standard_Integer ind = myDegree * aDimension; ind < (theDerivative + 1) * aDimension; ind++)
202 aPntDeriv[ind] = 0.0;
203 (&theDerivArray)[ind] = 0.0; // should be cleared separately, because aPntDeriv may look to another memory area
207 PLib::EvalPolynomial(aNewParameter, aDerivative, myDegree, aDimension,
208 aPolesArray[0], aPntDeriv[0]);
209 // Unnormalize derivatives since those are computed normalized
210 Standard_Real aFactor = 1.0;
211 for (Standard_Integer deriv = 1; deriv <= aDerivative; deriv++)
213 aFactor /= mySpanLength;
214 for (Standard_Integer ind = 0; ind < aDimension; ind++)
215 aPntDeriv[aDimension * deriv + ind] *= aFactor;
218 if (myIsRational) // calculate derivatives divided by weights derivatives
219 PLib::RationalDerivative(aDerivative, aDerivative, aDimension - 1, aPntDeriv[0], theDerivArray);
223 void BSplCLib_Cache::D0(const Standard_Real& theParameter, gp_Pnt2d& thePoint) const
225 Standard_Real aNewParameter = theParameter;
226 if (!myFlatKnots.IsNull()) // B-spline is periodical
227 PeriodicNormalization(myFlatKnots->Array1(), aNewParameter);
228 aNewParameter = (aNewParameter - mySpanStart) / mySpanLength;
230 Standard_Real* aPolesArray = ConvertArray(myPolesWeights);
231 Standard_Real aPoint[4];
232 Standard_Integer aDimension = myPolesWeights->RowLength(); // number of columns
234 PLib::NoDerivativeEvalPolynomial(aNewParameter, myDegree,
235 aDimension, myDegree * aDimension,
236 aPolesArray[0], aPoint[0]);
238 thePoint.SetCoord(aPoint[0], aPoint[1]);
240 thePoint.ChangeCoord().Divide(aPoint[2]);
243 void BSplCLib_Cache::D0(const Standard_Real& theParameter, gp_Pnt& thePoint) const
245 Standard_Real aNewParameter = theParameter;
246 if (!myFlatKnots.IsNull()) // B-spline is periodical
247 PeriodicNormalization(myFlatKnots->Array1(), aNewParameter);
248 aNewParameter = (aNewParameter - mySpanStart) / mySpanLength;
250 Standard_Real* aPolesArray = ConvertArray(myPolesWeights);
251 Standard_Real aPoint[4];
252 Standard_Integer aDimension = myPolesWeights->RowLength(); // number of columns
254 PLib::NoDerivativeEvalPolynomial(aNewParameter, myDegree,
255 aDimension, myDegree * aDimension,
256 aPolesArray[0], aPoint[0]);
258 thePoint.SetCoord(aPoint[0], aPoint[1], aPoint[2]);
260 thePoint.ChangeCoord().Divide(aPoint[3]);
264 void BSplCLib_Cache::D1(const Standard_Real& theParameter, gp_Pnt2d& thePoint, gp_Vec2d& theTangent) const
266 Standard_Integer aDimension = myPolesWeights->RowLength(); // number of columns
267 Standard_Real aPntDeriv[8]; // result storage (point and derivative coordinates)
269 this->CalculateDerivative(theParameter, 1, aPntDeriv[0]);
270 if (myIsRational) // the size of aPntDeriv was changed by PLib::RationalDerivative
273 thePoint.SetCoord(aPntDeriv[0], aPntDeriv[1]);
274 theTangent.SetCoord(aPntDeriv[aDimension], aPntDeriv[aDimension + 1]);
277 void BSplCLib_Cache::D1(const Standard_Real& theParameter, gp_Pnt& thePoint, gp_Vec& theTangent) const
279 Standard_Integer aDimension = myPolesWeights->RowLength(); // number of columns
280 Standard_Real aPntDeriv[8]; // result storage (point and derivative coordinates)
282 this->CalculateDerivative(theParameter, 1, aPntDeriv[0]);
283 if (myIsRational) // the size of aPntDeriv was changed by PLib::RationalDerivative
286 thePoint.SetCoord(aPntDeriv[0], aPntDeriv[1], aPntDeriv[2]);
287 theTangent.SetCoord(aPntDeriv[aDimension], aPntDeriv[aDimension + 1], aPntDeriv[aDimension + 2]);
290 void BSplCLib_Cache::D2(const Standard_Real& theParameter, gp_Pnt2d& thePoint, gp_Vec2d& theTangent, gp_Vec2d& theCurvature) const
292 Standard_Integer aDimension = myPolesWeights->RowLength(); // number of columns
293 Standard_Real aPntDeriv[12]; // result storage (point and derivatives coordinates)
295 this->CalculateDerivative(theParameter, 2, aPntDeriv[0]);
296 if (myIsRational) // the size of aPntDeriv was changed by PLib::RationalDerivative
299 thePoint.SetCoord(aPntDeriv[0], aPntDeriv[1]);
300 theTangent.SetCoord(aPntDeriv[aDimension], aPntDeriv[aDimension + 1]);
301 theCurvature.SetCoord(aPntDeriv[aDimension<<1], aPntDeriv[(aDimension<<1) + 1]);
304 void BSplCLib_Cache::D2(const Standard_Real& theParameter, gp_Pnt& thePoint, gp_Vec& theTangent, gp_Vec& theCurvature) const
306 Standard_Integer aDimension = myPolesWeights->RowLength(); // number of columns
307 Standard_Real aPntDeriv[12]; // result storage (point and derivatives coordinates)
309 this->CalculateDerivative(theParameter, 2, aPntDeriv[0]);
310 if (myIsRational) // the size of aPntDeriv was changed by PLib::RationalDerivative
313 thePoint.SetCoord(aPntDeriv[0], aPntDeriv[1], aPntDeriv[2]);
314 theTangent.SetCoord(aPntDeriv[aDimension], aPntDeriv[aDimension + 1], aPntDeriv[aDimension + 2]);
315 theCurvature.SetCoord(aPntDeriv[aDimension<<1], aPntDeriv[(aDimension<<1) + 1], aPntDeriv[(aDimension<<1) + 2]);
319 void BSplCLib_Cache::D3(const Standard_Real& theParameter,
321 gp_Vec2d& theTangent,
322 gp_Vec2d& theCurvature,
323 gp_Vec2d& theTorsion) const
325 Standard_Integer aDimension = myPolesWeights->RowLength(); // number of columns
326 Standard_Real aPntDeriv[16]; // result storage (point and derivatives coordinates)
328 this->CalculateDerivative(theParameter, 3, aPntDeriv[0]);
329 if (myIsRational) // the size of aPntDeriv was changed by PLib::RationalDerivative
332 thePoint.SetCoord(aPntDeriv[0], aPntDeriv[1]);
333 theTangent.SetCoord(aPntDeriv[aDimension], aPntDeriv[aDimension + 1]);
334 Standard_Integer aShift = aDimension << 1;
335 theCurvature.SetCoord(aPntDeriv[aShift], aPntDeriv[aShift + 1]);
336 aShift += aDimension;
337 theTorsion.SetCoord(aPntDeriv[aShift], aPntDeriv[aShift + 1]);
340 void BSplCLib_Cache::D3(const Standard_Real& theParameter,
343 gp_Vec& theCurvature,
344 gp_Vec& theTorsion) const
346 Standard_Integer aDimension = myPolesWeights->RowLength(); // number of columns
347 Standard_Real aPntDeriv[16]; // result storage (point and derivatives coordinates)
349 this->CalculateDerivative(theParameter, 3, aPntDeriv[0]);
350 if (myIsRational) // the size of aPntDeriv was changed by PLib::RationalDerivative
353 thePoint.SetCoord(aPntDeriv[0], aPntDeriv[1], aPntDeriv[2]);
354 theTangent.SetCoord(aPntDeriv[aDimension], aPntDeriv[aDimension + 1], aPntDeriv[aDimension + 2]);
355 Standard_Integer aShift = aDimension << 1;
356 theCurvature.SetCoord(aPntDeriv[aShift], aPntDeriv[aShift + 1], aPntDeriv[aShift + 2]);
357 aShift += aDimension;
358 theTorsion.SetCoord(aPntDeriv[aShift], aPntDeriv[aShift + 1], aPntDeriv[aShift + 2]);