projects / occt.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
0022048: Visualization, AIS_InteractiveContext - single object selection should alway...
[occt.git] / src / Extrema / Extrema_ExtPExtS.cxx
CommitLineData
b311480e 1// Created on: 1999-09-16
2// Created by: Edward AGAPOV
3// Copyright (c) 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.
7fd59977 16
17#include <Standard_NotImplemented.hxx>
18#include <Standard_OutOfRange.hxx>
19#include <StdFail_NotDone.hxx>
20#include <Adaptor3d_HCurve.hxx>
21#include <ElCLib.hxx>
22#include <Extrema_ExtPElC.hxx>
23#include <Extrema_ExtPExtS.hxx>
24#include <Extrema_POnCurv.hxx>
25#include <Extrema_POnSurf.hxx>
26#include <Precision.hxx>
27#include <gp.hxx>
28#include <gp_Ax2.hxx>
29#include <gp_Dir.hxx>
30#include <gp_Lin.hxx>
31#include <gp_Pln.hxx>
32#include <gp_Pnt.hxx>
33#include <gp_Vec.hxx>
34#include <math_FunctionSetRoot.hxx>
35#include <math_Vector.hxx>
6b84c3f7 36#include <GeomAdaptor_HSurfaceOfLinearExtrusion.hxx>
0734c53d 37
7fd59977 38
92efcf78 39IMPLEMENT_STANDARD_RTTIEXT(Extrema_ExtPExtS,Standard_Transient)
40
7fd59977 41static gp_Ax2 GetPosition (const Handle(Adaptor3d_HCurve)& C);
42
43static void PerformExtPElC (Extrema_ExtPElC& E,
44 const gp_Pnt& P,
45 const Handle(Adaptor3d_HCurve)& C,
46 const Standard_Real Tol);
47
48static Standard_Boolean
49 IsCaseAnalyticallyComputable (const GeomAbs_CurveType& theType,
50 const gp_Ax2& theCurvePos,
51 const gp_Dir& theSurfaceDirection) ;
52
53static gp_Pnt GetValue(const Standard_Real U,
54 const Handle(Adaptor3d_HCurve)& C);
55//=======================================================================
56//function : Project
57//purpose : Returns the projection of a point <Point> on a plane
58// <ThePlane> along a direction <TheDir>.
59//=======================================================================
60static gp_Pnt ProjectPnt(const gp_Ax2& ThePlane,
61 const gp_Dir& TheDir,
62 const gp_Pnt& Point)
63{
64 gp_Vec PO(Point,ThePlane.Location());
65 Standard_Real Alpha = PO * gp_Vec(ThePlane.Direction());
66 Alpha /= TheDir * ThePlane.Direction();
67 gp_Pnt P;
68 P.SetXYZ(Point.XYZ() + Alpha * TheDir.XYZ());
69 return P;
70}
71
72//=======================================================================
73//function : IsOriginalPnt
74//purpose :
75//=======================================================================
76
77static Standard_Boolean IsOriginalPnt (const gp_Pnt& P,
78 const Extrema_POnSurf* Points,
79 const Standard_Integer NbPoints)
80{
81 for (Standard_Integer i=1; i<=NbPoints; i++) {
82 if (Points[i-1].Value().IsEqual(P, Precision::Confusion())) {
83 return Standard_False;
84 }
85 }
86 return Standard_True;
87}
88
89//=======================================================================
90//function : MakePreciser
91//purpose :
92//=======================================================================
93
94void Extrema_ExtPExtS::MakePreciser (Standard_Real& U,
95 const gp_Pnt& P,
96 const Standard_Boolean isMin,
97 const gp_Ax2& OrtogSection) const
98{
99 if (U>myusup) {
100 U = myusup;
101 } else if (U<myuinf) {
102 U = myuinf;
103 } else {
104
105 Standard_Real step = (myusup - myuinf) / 30, D2e, D2next ,D2prev;
106 gp_Pnt
107 Pe = ProjectPnt (OrtogSection, myDirection, GetValue(U,myC)),
108 Pprev = ProjectPnt (OrtogSection, myDirection, GetValue(U-step, myC)),
109 Pnext = ProjectPnt (OrtogSection, myDirection, GetValue(U+step, myC));
110 D2e = P.SquareDistance(Pe),
111 D2next = P.SquareDistance(Pnext),
112 D2prev = P.SquareDistance(Pprev);
113 Standard_Boolean notFound;
114 if (isMin)
115 notFound = (D2e > D2prev || D2e > D2next);
116 else
117 notFound = (D2e < D2prev || D2e < D2next);
118
119 if (notFound && (D2e < D2next && isMin)) {
120 step = -step;
121 D2next = D2prev;
122 Pnext = Pprev;
123 }
124 while (notFound) {
125 U = U + step;
126 if (U > myusup) {
aabe3a17 127 U = myusup;
128 break;
7fd59977 129 }
130 if (U < myuinf) {
aabe3a17 131 U = myuinf;
132 break;
7fd59977 133 }
134 D2e = D2next;
135 Pe = Pnext;
136 Pnext = ProjectPnt (OrtogSection, myDirection, GetValue(U+step, myC));
137 D2next = P.SquareDistance(Pnext);
138 if (isMin)
aabe3a17 139 notFound = D2e > D2next;
7fd59977 140 else
aabe3a17 141 notFound = D2e < D2next;
7fd59977 142 }
143 }
144}
145//=============================================================================
146
cbff1e55 147Extrema_ExtPExtS::Extrema_ExtPExtS()
148: myuinf(0.0),
149 myusup(0.0),
150 mytolu(0.0),
151 myvinf(0.0),
152 myvsup(0.0),
153 mytolv(0.0),
154 myIsAnalyticallyComputable(Standard_False),
155 myDone(Standard_False),
156 myNbExt(Standard_False)
7fd59977 157{
7fd59977 158}
159
160//=============================================================================
161
0734c53d 162Extrema_ExtPExtS::Extrema_ExtPExtS (const gp_Pnt& theP,
6b84c3f7 163 const Handle(GeomAdaptor_HSurfaceOfLinearExtrusion)& theS,
0734c53d 164 const Standard_Real theUmin,
165 const Standard_Real theUsup,
166 const Standard_Real theVmin,
167 const Standard_Real theVsup,
168 const Standard_Real theTolU,
169 const Standard_Real theTolV)
cbff1e55 170: myuinf(theUmin),
171 myusup(theUsup),
172 mytolu(theTolU),
173 myvinf(theVmin),
174 myvsup(theVsup),
175 mytolv(theTolV),
176 myS (theS),
177 myIsAnalyticallyComputable(Standard_False),
178 myDone(Standard_False),
179 myNbExt(Standard_False)
7fd59977 180{
0734c53d 181 Initialize (theS,
182 theUmin,
183 theUsup,
184 theVmin,
185 theVsup,
186 theTolU,
187 theTolV);
188
189 Perform (theP);
7fd59977 190}
191//=============================================================================
192
0734c53d 193Extrema_ExtPExtS::Extrema_ExtPExtS (const gp_Pnt& theP,
6b84c3f7 194 const Handle(GeomAdaptor_HSurfaceOfLinearExtrusion)& theS,
0734c53d 195 const Standard_Real theTolU,
196 const Standard_Real theTolV)
cbff1e55 197: myuinf(theS->FirstUParameter()),
198 myusup(theS->LastUParameter()),
199 mytolu(theTolU),
200 myvinf(theS->FirstVParameter()),
201 myvsup(theS->LastVParameter()),
202 mytolv(theTolV),
203 myS (theS),
204 myIsAnalyticallyComputable(Standard_False),
205 myDone(Standard_False),
206 myNbExt(Standard_False)
7fd59977 207{
0734c53d 208 Initialize (theS,
209 theS->FirstUParameter(),
210 theS->LastUParameter(),
211 theS->FirstVParameter(),
212 theS->LastVParameter(),
213 theTolU,
214 theTolV);
215
216 Perform (theP);
7fd59977 217}
cbff1e55 218
7fd59977 219//=======================================================================
220//function : Initialize
221//purpose :
222//=======================================================================
223
6b84c3f7 224void Extrema_ExtPExtS::Initialize (const Handle(GeomAdaptor_HSurfaceOfLinearExtrusion)& theS,
0734c53d 225 const Standard_Real theUinf,
226 const Standard_Real theUsup,
227 const Standard_Real theVinf,
228 const Standard_Real theVsup,
229 const Standard_Real theTolU,
230 const Standard_Real theTolV)
7fd59977 231{
0734c53d 232 myuinf = theUinf;
233 myusup = theUsup;
234 mytolu = theTolU;
7fd59977 235
0734c53d 236 myvinf = theVinf;
237 myvsup = theVsup;
238 mytolv = theTolV;
7fd59977 239
0734c53d 240 Handle(Adaptor3d_HCurve) anACurve = theS->BasisCurve();
7fd59977 241
0734c53d 242 myF.Initialize (theS->ChangeSurface());
7fd59977 243 myC = anACurve;
0734c53d 244 myS = theS;
7fd59977 245 myPosition = GetPosition(myC);
0734c53d 246 myDirection = theS->Direction();
7fd59977 247 myIsAnalyticallyComputable = //Standard_False;
0734c53d 248 IsCaseAnalyticallyComputable (myC->GetType(), myPosition, myDirection);
7fd59977 249
250 if (!myIsAnalyticallyComputable)
0734c53d 251 {
252 myExtPS.Initialize (theS->ChangeSurface(),
253 32,
254 32,
255 theUinf,
256 theUsup,
257 theVinf,
258 theVsup,
259 theTolU,
260 theTolV);
261 }
7fd59977 262}
263
264
265//=======================================================================
266//function : Perform
267//purpose :
268//=======================================================================
269
270void Extrema_ExtPExtS::Perform (const gp_Pnt& P)
271{
aabe3a17 272 const Standard_Integer NbExtMax = 4; //dimension of arrays
273 //myPoint[] and mySqDist[]
274 //For "analytical" case
7fd59977 275 myDone = Standard_False;
276 myNbExt = 0;
277
278 if (!myIsAnalyticallyComputable) {
279 myExtPS.Perform(P);
280 myDone = myExtPS.IsDone();
281// modified by NIZHNY-EAP Wed Nov 17 12:59:08 1999 ___BEGIN___
282 myNbExt = myExtPS.NbExt();
283// modified by NIZHNY-EAP Wed Nov 17 12:59:09 1999 ___END___
284 return;
285 }
286
287 gp_Pnt Pe, Pp = ProjectPnt(myPosition,myDirection,P);
288 Extrema_ExtPElC anExt;
289 PerformExtPElC(anExt, Pp, myC, mytolu);
290 if (!anExt.IsDone()) return;
291
292 gp_Ax2 anOrtogSection (P, myDirection);
293 Standard_Real U,V;
294 Standard_Boolean
295 isMin,
296 isSimpleCase =
297 myDirection.IsParallel(myPosition.Direction(),Precision::Angular());
298 Standard_Integer i, aNbExt = anExt.NbExt();
299 math_Vector UV(1,2), Tol(1,2), UVinf(1,2), UVsup(1,2);
300 Tol(1) = mytolu; Tol(2) = mytolv;
301 UVinf(1) = myuinf; UVinf(2) = myvinf;
302 UVsup(1) = myusup; UVsup(2) = myvsup;
303
304
305 for (i=1; i<=aNbExt; i++) {
306 Extrema_POnCurv POC=anExt.Point(i);
307 U = POC.Parameter();
308 //// modified by jgv, 23.12.2008 for OCC17194 ////
309 if (myC->IsPeriodic())
aabe3a17 310 {
311 Standard_Real U2 = U;
312 ElCLib::AdjustPeriodic(myuinf, myuinf + 2.*M_PI, Precision::PConfusion(), U, U2);
313 }
7fd59977 314 //////////////////////////////////////////////////
315 gp_Pnt E = POC.Value();
316 Pe = ProjectPnt(anOrtogSection, myDirection, E);
317
318 if (isSimpleCase) {
319 V = gp_Vec(E,Pe) * gp_Vec(myDirection);
320 // modified by NIZHNY-MKK Thu Sep 18 14:46:14 2003.BEGIN
321 // myPoint[++myNbExt] = Extrema_POnSurf(U, V, Pe);
322 // myValue[myNbExt] = anExt.Value(i);
323 myPoint[myNbExt] = Extrema_POnSurf(U, V, Pe);
324 mySqDist[myNbExt] = anExt.SquareDistance(i);
325 myNbExt++;
aabe3a17 326 if(myNbExt == NbExtMax)
327 {
328 break;
329 }
7fd59977 330 // modified by NIZHNY-MKK Thu Sep 18 14:46:18 2003.END
331 }
332 else {
333 myF.SetPoint(P);
334 isMin = anExt.IsMin(i);//( Pp.Distance(GetValue(U+10,myC)) > anExt.Value(i) );
335
336 MakePreciser(U, P, isMin, anOrtogSection);
337 E = GetValue(U, myC);
338 Pe = ProjectPnt (anOrtogSection, myDirection, E),
339 V = gp_Vec(E,Pe) * gp_Vec(myDirection);
340 UV(1) = U; UV(2) = V;
859a47c3 341 math_FunctionSetRoot aFSR (myF, Tol);
342 aFSR.Perform(myF, UV, UVinf, UVsup);
7fd59977 343// for (Standard_Integer k=1 ; k <= myF.NbExt();
344 Standard_Integer k;
345 for ( k=1 ; k <= myF.NbExt(); k++) {
aabe3a17 346 if (IsOriginalPnt(myF.Point(k).Value(), myPoint, myNbExt)) {
347 // modified by NIZHNY-MKK Thu Sep 18 14:46:41 2003.BEGIN
348 // myPoint[++myNbExt] = myF.Point(k);
349 // myValue[myNbExt] = myF.Value(k);
350 myPoint[myNbExt] = myF.Point(k);
351 mySqDist[myNbExt] = myF.SquareDistance(k);
352 myNbExt++;
353 if(myNbExt == NbExtMax)
354 {
355 break;
356 }
357 // modified by NIZHNY-MKK Thu Sep 18 14:46:43 2003.END
358 }
359 }
360 if(myNbExt == NbExtMax)
361 {
362 break;
7fd59977 363 }
364 // try symmetric point
aabe3a17 365 myF.SetPoint(P); //To clear previous solutions
7fd59977 366 U *= -1;
367 MakePreciser(U, P, isMin, anOrtogSection);
368 E = GetValue(U, myC);
369 Pe = ProjectPnt (anOrtogSection, myDirection, E),
370 V = gp_Vec(E,Pe) * gp_Vec(myDirection);
371 UV(1) = U; UV(2) = V;
372
373 aFSR.Perform (myF,UV,UVinf,UVsup);
374
375 for (k=1 ; k <= myF.NbExt(); k++) {
aabe3a17 376 if(myF.SquareDistance(k) > Precision::Confusion()*Precision::Confusion())
377 {
378 //Additional checking solution: FSR sometimes is wrong
379 //when starting point is far from solution.
380 Standard_Real dist = Sqrt(myF.SquareDistance(k));
381 math_Vector Vals(1, 2);
382 const Extrema_POnSurf& PonS=myF.Point(k);
383 Standard_Real u, v;
384 PonS.Parameter(u, v);
385 UV(1) = u;
386 UV(2) = v;
387 myF.Value(UV, Vals);
388 gp_Vec du, dv;
389 myS->D1(u, v, Pe, du, dv);
390 Standard_Real mdu = du.Magnitude();
391 Standard_Real mdv = dv.Magnitude();
392 u = Abs(Vals(1));
393 v = Abs(Vals(2));
394 if(mdu > Precision::PConfusion())
395 {
396 if(u / dist / mdu > Precision::PConfusion())
397 {
398 continue;
399 }
400 }
401 if(mdv > Precision::PConfusion())
402 {
403 if(v / dist / mdv > Precision::PConfusion())
404 {
405 continue;
406 }
407 }
408
409 }
410 if (IsOriginalPnt(myF.Point(k).Value(), myPoint, myNbExt)) {
411 // modified by NIZHNY-MKK Thu Sep 18 14:46:59 2003.BEGIN
412 // myPoint[++myNbExt] = myF.Point(k);
413 // myValue[myNbExt] = myF.Value(k);
414 myPoint[myNbExt] = myF.Point(k);
415 mySqDist[myNbExt] = myF.SquareDistance(k);
416 myNbExt++;
417 if(myNbExt == NbExtMax)
418 {
419 break;
420 }
421 // modified by NIZHNY-MKK Thu Sep 18 14:47:04 2003.END
422 }
423 }
424 if(myNbExt == NbExtMax)
425 {
426 break;
7fd59977 427 }
428 }
429 }
430 myDone = Standard_True;
431 return;
432}
433
434//=============================================================================
435
436Standard_Boolean Extrema_ExtPExtS::IsDone () const { return myDone; }
437//=============================================================================
438
439Standard_Integer Extrema_ExtPExtS::NbExt () const
440{
9775fa61 441 if (!IsDone()) { throw StdFail_NotDone(); }
7fd59977 442 if (myIsAnalyticallyComputable)
443 return myNbExt;
444 else
445 return myExtPS.NbExt();
446}
447//=============================================================================
448
449Standard_Real Extrema_ExtPExtS::SquareDistance (const Standard_Integer N) const
450{
9775fa61 451 if (!IsDone()) { throw StdFail_NotDone(); }
452 if ((N < 1) || (N > myNbExt)) { throw Standard_OutOfRange(); }
7fd59977 453 if (myIsAnalyticallyComputable)
454 // modified by NIZHNY-MKK Thu Sep 18 14:48:39 2003.BEGIN
455 // return myValue[N];
456 return mySqDist[N-1];
457 // modified by NIZHNY-MKK Thu Sep 18 14:48:42 2003.END
458 else
459 return myExtPS.SquareDistance(N);
460}
461//=============================================================================
462
5d99f2c8 463const Extrema_POnSurf& Extrema_ExtPExtS::Point (const Standard_Integer N) const
7fd59977 464{
9775fa61 465 if (!IsDone()) { throw StdFail_NotDone(); }
466 if ((N < 1) || (N > myNbExt)) { throw Standard_OutOfRange(); }
7fd59977 467 if (myIsAnalyticallyComputable) {
468 // modified by NIZHNY-MKK Thu Sep 18 14:47:40 2003.BEGIN
469 // return myPoint[N];
470 return myPoint[N-1];
471 }
472 // modified by NIZHNY-MKK Thu Sep 18 14:47:43 2003.END
473 else
474 return myExtPS.Point(N);
475}
476//=============================================================================
477
478
479static gp_Ax2 GetPosition (const Handle(Adaptor3d_HCurve)& C)
480{
481 switch (C->GetType()) {
482 case GeomAbs_Line: {
483 gp_Lin L = C->Line();
484 gp_Pln Pln = gp_Pln(L.Location(), L.Direction());
485 //:abv 30.05.02: OCC - use constructor instead of Set...s() to avoid exception
486 gp_Ax2 Pos ( Pln.Location(), Pln.Position().Direction(), Pln.Position().XDirection() );
487// Pos.SetAxis(Pln.XAxis());
488// Pos.SetXDirection(Pln.YAxis().Direction());
489// Pos.SetYDirection(Pln.Position().Direction());
490 return Pos;
491 }
492 case GeomAbs_Circle:
493 return C->Circle().Position();
494 case GeomAbs_Ellipse:
495 return C->Ellipse().Position();
496 case GeomAbs_Hyperbola:
497 return C->Hyperbola().Position();
498 case GeomAbs_Parabola:
499 return C->Parabola().Position();
500 default:
501 return gp_Ax2 ();
502 }
503}
504//=============================================================================
505
506static void PerformExtPElC (Extrema_ExtPElC& E,
507 const gp_Pnt& P,
508 const Handle(Adaptor3d_HCurve)& C,
509 const Standard_Real Tol)
510{
511 switch (C->GetType()) {
512 case GeomAbs_Hyperbola:
513 E.Perform(P, C->Hyperbola(), Tol, -Precision::Infinite(),Precision::Infinite());
514 return;
515 case GeomAbs_Line:
516 E.Perform(P, C->Line(), Tol, -Precision::Infinite(),Precision::Infinite());
517 return;
518 case GeomAbs_Circle:
c6541a0c 519 E.Perform(P, C->Circle(), Tol, 0.0, 2.0 * M_PI);
7fd59977 520 return;
521 case GeomAbs_Ellipse:
c6541a0c 522 E.Perform(P, C->Ellipse(), Tol, 0.0, 2.0 * M_PI);
7fd59977 523 return;
524 case GeomAbs_Parabola:
525 E.Perform(P, C->Parabola(), Tol, -Precision::Infinite(),Precision::Infinite());
526 return;
7fd59977 527 default:
528 return;
7fd59977 529 }
530}
531
532//=======================================================================
533//function : IsCaseAnalyticallyComputable
534//purpose :
535//=======================================================================
536
537static Standard_Boolean IsCaseAnalyticallyComputable
538 (const GeomAbs_CurveType& theType,
539 const gp_Ax2& theCurvePos,
540 const gp_Dir& theSurfaceDirection)
541{
542 // check type
543 switch (theType) {
544 case GeomAbs_Line:
545 case GeomAbs_Circle:
546 case GeomAbs_Ellipse:
547 case GeomAbs_Hyperbola:
548 case GeomAbs_Parabola:
549 break;
550 default:
551 return Standard_False;
552 }
553 // check if it is a plane
554 if (Abs(theCurvePos.Direction() * theSurfaceDirection) <= gp::Resolution())
555 return Standard_False;
556 else
557 return Standard_True;
558}
559//=======================================================================
560//function : GetValue
561//purpose :
562//=======================================================================
563
564static gp_Pnt GetValue(const Standard_Real U,
565 const Handle(Adaptor3d_HCurve)& C)
566{
567 switch (C->GetType()) {
568 case GeomAbs_Line:
569 return ElCLib::Value(U, C->Line());
570 case GeomAbs_Circle:
571 return ElCLib::Value(U, C->Circle());
572 case GeomAbs_Ellipse:
573 return ElCLib::Value(U, C->Ellipse());
574 case GeomAbs_Hyperbola:
575 return ElCLib::Value(U, C->Hyperbola());
576 case GeomAbs_Parabola:
577 return ElCLib::Value(U, C->Parabola());
578 default:
579 return gp_Pnt();
580 }
581}
582//=======================================================================
583//function : GetU
584//purpose :
585//=======================================================================
0797d9d3 586//#ifdef OCCT_DEBUG
7fd59977 587//static Standard_Real GetU(const gp_Vec& vec,
588// const gp_Pnt& P,
589// const Handle(Adaptor3d_HCurve)& C)
590//{
591// switch (C->GetType()) {
592// case GeomAbs_Line:
593// return ElCLib::Parameter(C->Line().Translated(vec), P);
594// case GeomAbs_Circle:
595// return ElCLib::Parameter(C->Circle().Translated(vec), P);
596// case GeomAbs_Ellipse:
597// return ElCLib::Parameter(C->Ellipse().Translated(vec), P);
598// case GeomAbs_Hyperbola:
599// return ElCLib::Parameter(C->Hyperbola().Translated(vec), P);
600// case GeomAbs_Parabola:
601// return ElCLib::Parameter(C->Parabola().Translated(vec), P);
602// default:
603// return 0;
604// }
605//}
606//#endif
Open CASCADE Technology repository