0026605: Possible array out of bounds read in Extrema_GExtPC.gxx
[occt.git] / src / Extrema / Extrema_GExtPC.gxx
CommitLineData
b311480e 1// Created on: 1992-10-19
2// Created by: Laurent PAINNOT
3// Copyright (c) 1992-1999 Matra Datavision
973c2be1 4// Copyright (c) 1999-2014 OPEN CASCADE SAS
b311480e 5//
973c2be1 6// This file is part of Open CASCADE Technology software library.
b311480e 7//
d5f74e42 8// This library is free software; you can redistribute it and/or modify it under
9// the terms of the GNU Lesser General Public License version 2.1 as published
973c2be1 10// by the Free Software Foundation, with special exception defined in the file
11// OCCT_LGPL_EXCEPTION.txt. Consult the file LICENSE_LGPL_21.txt included in OCCT
12// distribution for complete text of the license and disclaimer of any warranty.
b311480e 13//
973c2be1 14// Alternatively, this file may be used under the terms of Open CASCADE
15// commercial license or contractual agreement.
b311480e 16
7fd59977 17// Modified by cma, Fri Dec 10 18:11:56 1993
7fd59977 18
19#include Extrema_EPC_hxx
20#include ThePOnC_hxx
21#include TheSequenceOfPOnC_hxx
22#include ThePoint_hxx
23#include TheVector_hxx
24#include TheExtPElC_hxx
25#include <StdFail_NotDone.hxx>
26#include <Standard_Failure.hxx>
27#include <GeomAbs_CurveType.hxx>
28#include <Precision.hxx>
29#include <ElCLib.hxx>
30#include <TColStd_Array1OfReal.hxx>
c8bf1eb7 31#include <NCollection_Array1.hxx>
7fd59977 32
33
34//=======================================================================
35//function : Perform
36//purpose :
37//=======================================================================
7fd59977 38void Extrema_GExtPC::Perform(const ThePoint& P)
39{
40 mySqDist.Clear();
41 mypoint.Clear();
42 myismin.Clear();
43 Standard_Integer i, NbExt, n;
44 Standard_Real U;
45 mysample = 17;
46 Standard_Real t3d = Precision::Confusion();
47 if (Precision::IsInfinite(myuinf)) mydist1 = RealLast();
48 else {
49 Pf = TheCurveTool::Value(*((TheCurve*)myC), myuinf);
50 mydist1 = P.SquareDistance(Pf);
51 }
52
53 if (Precision::IsInfinite(myusup)) mydist2 = RealLast();
54 else {
55 Pl = TheCurveTool::Value(*((TheCurve*)myC), myusup);
56 mydist2 = P.SquareDistance(Pl);
57 }
58
c8bf1eb7 59 TheCurve & aCurve = *((TheCurve*)myC);
60
7fd59977 61 switch(type) {
62 case GeomAbs_Circle:
63 {
c8bf1eb7 64 myExtPElC.Perform(P, TheCurveTool::Circle(aCurve), t3d, myuinf, myusup);
7fd59977 65 }
66 break;
67 case GeomAbs_Ellipse:
68 {
c8bf1eb7 69 myExtPElC.Perform(P, TheCurveTool::Ellipse(aCurve), t3d, myuinf, myusup);
7fd59977 70 }
71 break;
72 case GeomAbs_Parabola:
73 {
c8bf1eb7 74 myExtPElC.Perform(P, TheCurveTool::Parabola(aCurve), t3d,myuinf,myusup);
7fd59977 75 }
76 break;
77 case GeomAbs_Hyperbola:
78 {
c8bf1eb7 79 myExtPElC.Perform(P,TheCurveTool::Hyperbola(aCurve),t3d, myuinf, myusup);
7fd59977 80 }
81 break;
82 case GeomAbs_Line:
83 {
c8bf1eb7 84 myExtPElC.Perform(P, TheCurveTool::Line(aCurve), t3d, myuinf, myusup);
7fd59977 85 }
86 break;
87 case GeomAbs_BezierCurve:
88 {
89 myintuinf = myuinf;
90 myintusup = myusup;
c8bf1eb7 91 mysample = (TheCurveTool::Bezier(aCurve))->NbPoles() * 2;
92 myExtPC.Initialize(aCurve);
7fd59977 93 IntervalPerform(P);
94 return;
95 }
96 case GeomAbs_BSplineCurve:
97 {
c8bf1eb7 98 const Standard_Integer
99 aFirstIdx = TheCurveTool::BSpline(aCurve)->FirstUKnotIndex(),
100 aLastIdx = TheCurveTool::BSpline(aCurve)->LastUKnotIndex();
101 // const reference can not be used due to implementation of BRep_Adaptor.
102 TColStd_Array1OfReal aKnots(aFirstIdx, aLastIdx);
103 TheCurveTool::BSpline(aCurve)->Knots(aKnots);
104
105 // Workaround to work with:
106 // blend, where knots may be moved from param space.
107 Standard_Real aPeriodJump = 0.0;
0cbfb9f1 108 // Avoid prolem with too close knots.
109 const Standard_Real aTolCoeff = (myusup - myuinf) * Precision::PConfusion();
c8bf1eb7 110 if (TheCurveTool::IsPeriodic(aCurve))
111 {
112 Standard_Integer aPeriodShift =
113 Standard_Integer ((myuinf - aKnots(aFirstIdx)) / TheCurveTool::Period(aCurve));
114
0cbfb9f1 115 if (myuinf < aKnots(aFirstIdx) - aTolCoeff)
c8bf1eb7 116 aPeriodShift--;
117
118 aPeriodJump = TheCurveTool::Period(aCurve) * aPeriodShift;
119 }
120
121 Standard_Integer anIdx;
122
1581d651 123 // Find first and last used knot.
c8bf1eb7 124 Standard_Integer aFirstUsedKnot = aFirstIdx,
125 aLastUsedKnot = aLastIdx;
126 for(anIdx = aFirstIdx; anIdx <= aLastIdx; anIdx++)
127 {
128 Standard_Real aKnot = aKnots(anIdx) + aPeriodJump;
0cbfb9f1 129 if (myuinf >= aKnot - aTolCoeff)
c8bf1eb7 130 aFirstUsedKnot = anIdx;
131 else
132 break;
133
134 }
135 for(anIdx = aLastIdx; anIdx >= aFirstIdx; anIdx--)
136 {
137 Standard_Real aKnot = aKnots(anIdx) + aPeriodJump;
0cbfb9f1 138 if (myusup <= aKnot + aTolCoeff)
c8bf1eb7 139 aLastUsedKnot = anIdx;
140 else
141 break;
142 }
143
1581d651 144 if (aFirstUsedKnot == aLastUsedKnot)
145 {
146 // Degenerated case:
147 // Some bounds lies out of curve param space.
148 // In this case build one interval with [myuinf, myusup].
149 // Parameters of these indexes will be redefined.
150 aFirstUsedKnot = aFirstIdx;
151 aLastUsedKnot = aFirstIdx + 1;
152 }
153
c8bf1eb7 154 mysample = (TheCurveTool::BSpline(aCurve))->Degree() + 1;
155
156 // Fill sample points.
157 Standard_Integer aValIdx = 1;
158 NCollection_Array1<Standard_Real> aVal(1, (mysample) * (aLastUsedKnot - aFirstUsedKnot) + 1);
159 NCollection_Array1<Standard_Real> aParam(1, (mysample) * (aLastUsedKnot - aFirstUsedKnot) + 1);
160 for(anIdx = aFirstUsedKnot; anIdx < aLastUsedKnot; anIdx++)
161 {
162 Standard_Real aF = aKnots(anIdx) + aPeriodJump,
163 aL = aKnots(anIdx + 1) + aPeriodJump;
164
165 if (anIdx == aFirstUsedKnot)
166 aF = myuinf;
167 if (anIdx == aLastUsedKnot - 1)
168 aL = myusup;
169
170 Standard_Real aStep = (aL - aF) / mysample;
171 for(Standard_Integer aPntIdx = 0; aPntIdx < mysample; aPntIdx++)
172 {
173 Standard_Real aCurrentParam = aF + aStep * aPntIdx;
174 aVal(aValIdx) = TheCurveTool::Value(aCurve, aCurrentParam).SquareDistance(P);
175 aParam(aValIdx) = aCurrentParam;
176 aValIdx++;
177 }
178 }
179 // Fill last point.
180 aVal(aValIdx) = TheCurveTool::Value(aCurve, myusup).SquareDistance(P);
181 aParam(aValIdx) = myusup;
182
183 myExtPC.Initialize(aCurve);
184
185 // Find extremas.
186 for(anIdx = aVal.Lower() + 1; anIdx < aVal.Upper(); anIdx++)
187 {
188 if (aVal(anIdx) <= Precision::SquareConfusion())
189 {
190 mySqDist.Append(aVal(anIdx));
191 myismin.Append(Standard_True);
192 mypoint.Append(ThePOnC(aParam(anIdx), TheCurveTool::Value(aCurve, aParam(anIdx))));
193 }
194 if ((aVal(anIdx) >= aVal(anIdx + 1) &&
195 aVal(anIdx) >= aVal(anIdx - 1)) ||
196 (aVal(anIdx) <= aVal(anIdx + 1) &&
197 aVal(anIdx) <= aVal(anIdx - 1)) )
198 {
199 myintuinf = aParam(anIdx - 1);
200 myintusup = aParam(anIdx + 1);
201
202 IntervalPerform(P);
203 }
204 }
205
206 // Solve on first and last interval.
207 if (mydist1 > Precision::SquareConfusion())
208 {
209 ThePoint aP1, aP2;
210 TheVector aV1, aV2;
211 TheCurveTool::D1(aCurve, aParam.Value(aParam.Lower()), aP1, aV1);
212 TheCurveTool::D1(aCurve, aParam.Value(aParam.Lower() + 1), aP2, aV2);
213 TheVector aBase1(P, aP1), aBase2(P, aP2);
214 Standard_Real aVal1 = aV1.Dot(aBase1); // Derivative of (C(u) - P)^2
215 Standard_Real aVal2 = aV2.Dot(aBase2); // Derivative of (C(u) - P)^2
216
217 // Derivatives have opposite signs - min or max inside of interval (sufficient condition).
218 // Necessary condition - when point lies on curve.
219 if(aVal1 * aVal2 <= 0.0 ||
220 aBase1.Dot(aBase2) <= 0.0)
221 {
222 myintuinf = aParam(aVal.Lower());
223 myintusup = aParam(aVal.Lower() + 1);
224 IntervalPerform(P);
225 }
226 }
227
228 if (mydist2 > Precision::SquareConfusion())
229 {
230 ThePoint aP1, aP2;
231 TheVector aV1, aV2;
232 TheCurveTool::D1(aCurve, aParam.Value(aParam.Upper() - 1), aP1, aV1);
233 TheCurveTool::D1(aCurve, aParam.Value(aParam.Upper()), aP2, aV2);
234 TheVector aBase1(P, aP1), aBase2(P, aP2);
235 Standard_Real aVal1 = aV1.Dot(aBase1); // Derivative of (C(u) - P)^2
236 Standard_Real aVal2 = aV2.Dot(aBase2); // Derivative of (C(u) - P)^2
237
238 // Derivatives have opposite signs - min or max inside of interval (sufficient condition).
239 // Necessary condition - when point lies on curve.
240 if(aVal1 * aVal2 <= 0.0 ||
241 aBase1.Dot(aBase2) <= 0.0)
242 {
243 myintuinf = aParam(aVal.Upper() - 1);
244 myintusup = aParam(aVal.Upper());
245 IntervalPerform(P);
246 }
247 }
248
249 mydone = Standard_True;
250 break;
7fd59977 251 }
c8bf1eb7 252 case GeomAbs_OtherCurve:
7fd59977 253 {
254 Standard_Boolean IntExtIsDone = Standard_False;
255 Standard_Boolean IntIsNotValid;
c8bf1eb7 256 n = TheCurveTool::NbIntervals(aCurve, GeomAbs_C2);
7fd59977 257 TColStd_Array1OfReal theInter(1, n+1);
c8bf1eb7 258 Standard_Boolean isPeriodic = TheCurveTool::IsPeriodic(aCurve);
259 TheCurveTool::Intervals(aCurve, theInter, GeomAbs_C2);
7fd59977 260 mysample = Max(mysample/n, 17);
261 TheVector V1;
262 ThePoint PP;
263 Standard_Real s1 = 0.0 ;
c8bf1eb7 264 Standard_Real s2 = 0.0;
265 myExtPC.Initialize(aCurve);
266 for (i = 1; i <= n; i++)
267 {
253881cf 268 myintuinf = theInter(i);
269 myintusup = theInter(i+1);
270
271 Standard_Real anInfToCheck = myintuinf;
272 Standard_Real aSupToCheck = myintusup;
273
274 if (isPeriodic) {
c8bf1eb7 275 Standard_Real aPeriod = TheCurveTool::Period(aCurve);
253881cf 276 anInfToCheck = ElCLib::InPeriod(myintuinf, myuinf, myuinf+aPeriod);
277 aSupToCheck = myintusup+(anInfToCheck-myintuinf);
278 }
279 IntIsNotValid = (myuinf > aSupToCheck) || (myusup < anInfToCheck);
280
281 if(IntIsNotValid) continue;
282
283 if (myuinf >= anInfToCheck) anInfToCheck = myuinf;
284 if (myusup <= aSupToCheck) aSupToCheck = myusup;
285 if((aSupToCheck - anInfToCheck) <= mytolu) continue;
32ca7a51 286
287 if (i != 1)
288 {
c8bf1eb7 289 TheCurveTool::D1(aCurve, myintuinf, PP, V1);
253881cf 290 s1 = (TheVector(P, PP))*V1;
291 if (s1*s2 < 0.0) {
292 mySqDist.Append(PP.SquareDistance(P));
293 myismin.Append((s1 < 0.0));
294 mypoint.Append(ThePOnC(myintuinf, PP));
295 }
296 }
297 if (i != n) {
c8bf1eb7 298 TheCurveTool::D1(aCurve, myintusup, PP, V1);
253881cf 299 s2 = (TheVector(P, PP))*V1;
300 }
32ca7a51 301
253881cf 302 IntervalPerform(P);
303 IntExtIsDone = IntExtIsDone || mydone;
7fd59977 304 }
c8bf1eb7 305
7fd59977 306 mydone = IntExtIsDone;
c8bf1eb7 307 break;
308 }
309 }
94f71cad 310
c8bf1eb7 311 // Postprocessing.
312 if (type == GeomAbs_BSplineCurve ||
313 type == GeomAbs_OtherCurve)
314 {
315 // Additional checking if the point is on the first or last point of the curve
316 // and does not added yet.
317 if (mydist1 < Precision::SquareConfusion() ||
318 mydist2 < Precision::SquareConfusion())
319 {
320 Standard_Boolean isFirstAdded = Standard_False;
321 Standard_Boolean isLastAdded = Standard_False;
322 Standard_Integer aNbPoints = mypoint.Length();
323 for (i = 1; i <= aNbPoints; i++)
94f71cad 324 {
c8bf1eb7 325 U = mypoint.Value(i).Parameter();
326 if (Abs(U - myuinf) < mytolu)
327 isFirstAdded = Standard_True;
328 else if (Abs(myusup - U) < mytolu)
329 isLastAdded = Standard_True;
330 }
331 if (!isFirstAdded && mydist1 < Precision::SquareConfusion())
332 {
333 mySqDist.Prepend(mydist1);
334 myismin.Prepend(Standard_True);
335 mypoint.Prepend(ThePOnC(myuinf, Pf));
336 }
337 if (!isLastAdded && mydist2 < Precision::SquareConfusion())
338 {
339 mySqDist.Append(mydist2);
340 myismin.Append(Standard_True);
341 mypoint.Append(ThePOnC(myusup, Pl));
342 }
343 mydone = Standard_True;
344 }
345 }
346 else
347 {
348 // In analytical case
349 mydone = myExtPElC.IsDone();
350 if (mydone)
351 {
352 NbExt = myExtPElC.NbExt();
353 for (i = 1; i <= NbExt; i++)
354 {
355 // Verification de la validite des parametres:
356 ThePOnC PC = myExtPElC.Point(i);
357 U = PC.Parameter();
358 if (TheCurveTool::IsPeriodic(aCurve))
94f71cad 359 {
c8bf1eb7 360 U = ElCLib::InPeriod(U, myuinf, myuinf+TheCurveTool::Period(aCurve));
94f71cad 361 }
c8bf1eb7 362 if ((U >= myuinf-mytolu) && (U <= myusup+mytolu))
94f71cad 363 {
c8bf1eb7 364 PC.SetValues(U, myExtPElC.Point(i).Value());
365 mySqDist.Append(myExtPElC.SquareDistance(i));
366 myismin.Append(myExtPElC.IsMin(i));
367 mypoint.Append(PC);
94f71cad 368 }
369 }
7fd59977 370 }
371 }
7fd59977 372}
373
374
375//=======================================================================
376//function : Initialize
377//purpose :
378//=======================================================================
379
380void Extrema_GExtPC::Initialize(const TheCurve& C,
381 const Standard_Real Uinf,
382 const Standard_Real Usup,
383 const Standard_Real TolF)
384{
385 myC = (Standard_Address)&C;
386 myintuinf = myuinf = Uinf;
387 myintusup = myusup = Usup;
388 mytolf = TolF;
389 mytolu = TheCurveTool::Resolution(*((TheCurve*)myC), Precision::Confusion());
390 type = TheCurveTool::GetType(C);
391 mydone = Standard_False;
392 mydist1 = RealLast();
393 mydist2 = RealLast();
394 mysample = 17;
395}
396
397
398//=======================================================================
399//function : IntervalPerform
400//purpose :
401//=======================================================================
402
403void Extrema_GExtPC::IntervalPerform(const ThePoint& P)
404{
405 Standard_Integer i;
406 Standard_Real U;
c8bf1eb7 407 myExtPC.Initialize(mysample, myintuinf, myintusup, mytolu, mytolf);
7fd59977 408 myExtPC.Perform(P);
409 mydone = myExtPC.IsDone();
c8bf1eb7 410 if (mydone)
411 {
7fd59977 412 Standard_Integer NbExt = myExtPC.NbExt();
c8bf1eb7 413 for (i = 1; i <= NbExt; i++)
414 {
7fd59977 415 // Verification de la validite des parametres pour le cas trimme:
416 ThePOnC PC = myExtPC.Point(i);
417 U = PC.Parameter();
c8bf1eb7 418 if (TheCurveTool::IsPeriodic(*((TheCurve*)myC)))
419 {
253881cf 420 U = ElCLib::InPeriod(U, myuinf, myuinf+TheCurveTool::Period(*((TheCurve*)myC)));
7fd59977 421 }
c8bf1eb7 422 if ((U >= myuinf - mytolu) && (U <= myusup + mytolu))
423 {
253881cf 424 PC.SetValues(U, PC.Value());
425 mySqDist.Append(myExtPC.SquareDistance(i));
426 myismin.Append(myExtPC.IsMin(i));
427 mypoint.Append(PC);
7fd59977 428 }
429 }
430 }
431}
432
433
434
435
436//=======================================================================
437//function : Extrema_GExtPC
438//purpose :
439//=======================================================================
440
441Extrema_GExtPC::Extrema_GExtPC()
442{
443 myC = 0;
444 mydone = Standard_False;
445 mydist1 = RealLast();
446 mydist2 = RealLast();
447 mytolu = 0.0;
448 mytolf = 0.0;
449 mysample = 17;
450 myintuinf = myintusup = myuinf = myusup = Precision::Infinite();
451 type = GeomAbs_OtherCurve;
452}
453
454//=======================================================================
455//function : Extrema_GExtPC
456//purpose :
457//=======================================================================
458
459Extrema_GExtPC::Extrema_GExtPC(const ThePoint& P,
460 const TheCurve& C,
461 const Standard_Real Uinf,
462 const Standard_Real Usup,
463 const Standard_Real TolF)
464{
465 Initialize(C, Uinf, Usup, TolF);
466 Perform(P);
467}
468
469//=======================================================================
470//function : Extrema_GExtPC
471//purpose :
472//=======================================================================
473
474Extrema_GExtPC::Extrema_GExtPC(const ThePoint& P,
475 const TheCurve& C,
476 const Standard_Real TolF)
477{
478 Initialize(C, TheCurveTool::FirstParameter(C),
479 TheCurveTool::LastParameter(C), TolF);
480 Perform(P);
481}
482
483
484//=======================================================================
485//function : IsDone
486//purpose :
487//=======================================================================
488
489Standard_Boolean Extrema_GExtPC::IsDone() const
490{
491 return mydone;
492}
493
494
495//=======================================================================
496//function : Value
497//purpose :
498//=======================================================================
499
500Standard_Real Extrema_GExtPC::SquareDistance(const Standard_Integer N) const
501{
502 if(!mydone) StdFail_NotDone::Raise();
503 if ((N < 1) || (N > mySqDist.Length())) Standard_OutOfRange::Raise();
504 return mySqDist.Value(N);
505}
506
507
508//=======================================================================
509//function : NbExt
510//purpose :
511//=======================================================================
512
513Standard_Integer Extrema_GExtPC::NbExt() const
514{
515 if(!mydone) StdFail_NotDone::Raise();
516 return mySqDist.Length();
517}
518
519
520//=======================================================================
521//function : IsMin
522//purpose :
523//=======================================================================
524
525Standard_Boolean Extrema_GExtPC::IsMin(const Standard_Integer N) const
526{
527 if(!mydone) StdFail_NotDone::Raise();
528 if ((N < 1) || (N > mySqDist.Length())) Standard_OutOfRange::Raise();
529 return myismin.Value(N);
530}
531
532
533
534//=======================================================================
535//function : Point
536//purpose :
537//=======================================================================
538
5d99f2c8 539const ThePOnC & Extrema_GExtPC::Point(const Standard_Integer N) const
7fd59977 540{
541 if(!mydone) StdFail_NotDone::Raise();
542 if ((N < 1) || (N > mySqDist.Length())) Standard_OutOfRange::Raise();
543 return mypoint.Value(N);
544}
545
546
547//=======================================================================
548//function : TrimmedDistances
549//purpose :
550//=======================================================================
551
552void Extrema_GExtPC::TrimmedSquareDistances(Standard_Real& dist1,
553 Standard_Real& dist2,
554 ThePoint& P1,
555 ThePoint& P2) const
556{
557 dist1 = mydist1;
558 dist2 = mydist2;
559 P1 = Pf;
560 P2 = Pl;
561}