projects / occt.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
0028966: Coding Rules - remove Adaptor2d_HCurve2d, Adaptor3d_HCurve and Adaptor3d_HSu...
[occt.git] / src / IntStart / IntStart_SearchOnBoundaries.gxx
CommitLineData
b311480e 1// Copyright (c) 1995-1999 Matra Datavision
973c2be1 2// Copyright (c) 1999-2014 OPEN CASCADE SAS
b311480e 3//
973c2be1 4// This file is part of Open CASCADE Technology software library.
b311480e 5//
d5f74e42 6// This library is free software; you can redistribute it and/or modify it under
7// the terms of the GNU Lesser General Public License version 2.1 as published
973c2be1 8// by the Free Software Foundation, with special exception defined in the file
9// OCCT_LGPL_EXCEPTION.txt. Consult the file LICENSE_LGPL_21.txt included in OCCT
10// distribution for complete text of the license and disclaimer of any warranty.
b311480e 11//
973c2be1 12// Alternatively, this file may be used under the terms of Open CASCADE
13// commercial license or contractual agreement.
7fd59977 14
f542b7bb 15#include <algorithm>
9369e98a 16#include <memory>
74d80fb9 17#include <TopoDS_Edge.hxx>
18#include <Geom_Curve.hxx>
19#include <BRepAdaptor_Curve.hxx>
c22b52d6 20#include <Adaptor3d_Surface.hxx>
21#include <Adaptor3d_CurveOnSurface.hxx>
9369e98a 22#include <Adaptor3d_CurveOnSurface.hxx>
74d80fb9 23#include <GeomAbs_SurfaceType.hxx>
24#include <BRep_Tool.hxx>
25#include <Geom_Line.hxx>
9369e98a 26#include <Geom_Plane.hxx>
27#include <Geom_CylindricalSurface.hxx>
28#include <Geom_ConicalSurface.hxx>
29#include <Geom_SphericalSurface.hxx>
30#include <Geom_ToroidalSurface.hxx>
74d80fb9 31#include <gp_Lin.hxx>
32#include <gp_Vec.hxx>
33#include <gp_Dir.hxx>
34#include <gp_Cylinder.hxx>
35#include <gp_Ax1.hxx>
36#include <gp_Lin.hxx>
37
38#include <GeomAdaptor_Curve.hxx>
c22b52d6 39#include <GeomAdaptor_Surface.hxx>
74d80fb9 40#include <Precision.hxx>
41#include <Extrema_ExtCC.hxx>
9369e98a 42//#include <Extrema_ExtCS.hxx>
74d80fb9 43#include <Extrema_POnCurv.hxx>
9369e98a 44#include <IntCurveSurface_HInter.hxx>
74d80fb9 45
7fd59977 46#include <math_FunctionSample.hxx>
47#include <math_FunctionAllRoots.hxx>
7fd59977 48#include <TColgp_SequenceOfPnt.hxx>
49
74d80fb9 50// Modified by skv - Tue Aug 31 12:13:51 2004 OCC569
51
52#include <Precision.hxx>
53#include <IntSurf_Quadric.hxx>
f542b7bb 54#include <math_Function.hxx>
55#include <math_BrentMinimum.hxx>
56#include <math_Matrix.hxx>
57#include <math_Vector.hxx>
58#include <NCollection_Array1.hxx>
74d80fb9 59
9369e98a 60#ifdef OCCT_DEBUG
61#include <Geom_Circle.hxx>
62#include <Geom_Ellipse.hxx>
63#include <Geom_Hyperbola.hxx>
64#include <Geom_Parabola.hxx>
65#include <Geom_BezierCurve.hxx>
66#include <Geom_BSplineCurve.hxx>
67#include <GeomLib.hxx>
68#endif
69
70
c22b52d6 71static Standard_Boolean IsDegenerated(const Handle(Adaptor3d_CurveOnSurface)& theCurve);
9369e98a 72static Standard_Boolean IsDegenerated(const IntSurf_Quadric& theQuadric);
73
74d80fb9 74static void FindVertex (const TheArc&,
75 const Handle(TheTopolTool)&,
76 TheFunction&,
77 IntStart_SequenceOfPathPoint&,
78 const Standard_Real);
79
80
81static void BoundedArc (const TheArc& A,
82 const Handle(TheTopolTool)& Domain,
83 const Standard_Real Pdeb,
84 const Standard_Real Pfin,
85 TheFunction& Func,
86 IntStart_SequenceOfPathPoint& pnt,
87 IntStart_SequenceOfSegment& seg,
88 const Standard_Real TolBoundary,
89 const Standard_Real TolTangency,
90 Standard_Boolean& Arcsol,
91 const Standard_Boolean RecheckOnRegularity);
92
74d80fb9 93static void PointProcess (const gp_Pnt&,
94 const Standard_Real,
95 const TheArc&,
96 const Handle(TheTopolTool)&,
97 IntStart_SequenceOfPathPoint&,
98 const Standard_Real,
99 Standard_Integer&);
100
101static Standard_Integer TreatLC (const TheArc& A,
102 const Handle(TheTopolTool)& aDomain,
103 const IntSurf_Quadric& aQuadric,
104 const Standard_Real TolBoundary,
105 IntStart_SequenceOfPathPoint& pnt);
106
107static Standard_Boolean IsRegularity(const TheArc& A,
108 const Handle(TheTopolTool)& aDomain);
109
f542b7bb 110class MinFunction : public math_Function
111{
112public:
113 MinFunction(TheFunction &theFunc) : myFunc(&theFunc) {};
114
115 //returns value of the one-dimension-function when parameter
116 //is equal to theX
117 virtual Standard_Boolean Value(const Standard_Real theX,
118 Standard_Real& theFVal)
119 {
120 if(!myFunc->Value(theX, theFVal))
121 return Standard_False;
122
123 theFVal *= theFVal;
124 return Standard_True;
125 }
126
127 //see analogical method for abstract owner class math_Function
128 virtual Standard_Integer GetStateNumber()
129 {
130 return 0;
131 }
132
133private:
134 TheFunction *myFunc;
135};
136
74d80fb9 137
138//=======================================================================
139//function : FindVertex
140//purpose :
141//=======================================================================
142void FindVertex (const TheArc& A,
143 const Handle(TheTopolTool)& Domain,
144 TheFunction& Func,
145 IntStart_SequenceOfPathPoint& pnt,
146 const Standard_Real Toler)
147{
148
149// Find the vertex of the arc A restriction solutions. It stores
150// Vertex in the list solutions pnt.
151
152
153 TheVertex vtx;
154 Standard_Real param,valf;
155 Standard_Integer itemp;
156
157 Domain->Initialize(A);
158 Domain->InitVertexIterator();
159 while (Domain->MoreVertex()) {
160 vtx = Domain->Vertex();
161 param = TheSOBTool::Parameter(vtx,A);
162
163 // Evaluate the function and look compared to tolerance of the
164 // Vertex. If distance <= tolerance then add a vertex to the list of solutions.
165 // The arc is already assumed in the load function.
166
167 Func.Value(param,valf);
168 if (Abs(valf) <= Toler) {
169 itemp = Func.GetStateNumber();
170 pnt.Append(IntStart_ThePathPoint(Func.Valpoint(itemp),Toler, vtx,A,param));
171 // Solution is added
172 }
173 Domain->NextVertex();
174 }
175}
176
c22b52d6 177Standard_Boolean IsDegenerated(const Handle(Adaptor3d_CurveOnSurface)& theCurve)
9369e98a 178{
179 if (theCurve->GetType() == GeomAbs_Circle)
180 {
181 gp_Circ aCirc = theCurve->Circle();
182 if (aCirc.Radius() <= Precision::Confusion())
183 return Standard_True;
184 }
185 return Standard_False;
186}
187
188Standard_Boolean IsDegenerated(const IntSurf_Quadric& theQuadric)
189{
190 GeomAbs_SurfaceType TypeQuad = theQuadric.TypeQuadric();
191 if (TypeQuad == GeomAbs_Cone)
192 {
193 gp_Cone aCone = theQuadric.Cone();
194 Standard_Real aSemiAngle = Abs(aCone.SemiAngle());
195 if (aSemiAngle < 0.02 || aSemiAngle > 1.55)
196 return Standard_True;
197 }
198 return Standard_False;
199}
200
f542b7bb 201class SolInfo
202{
203public:
204 SolInfo() : myMathIndex(-1), myValue(RealLast())
205 {
206 }
207
208 void Init(const math_FunctionAllRoots& theSolution, const Standard_Integer theIndex)
209 {
210 myMathIndex = theIndex;
211 myValue = theSolution.GetPoint(theIndex);
212 }
213
9369e98a 214 void Init(const IntCurveSurface_HInter& theSolution, const Standard_Integer theIndex)
215 {
216 myMathIndex = theIndex;
217 myValue = theSolution.Point(theIndex).W();
218 }
219
f542b7bb 220 Standard_Real Value() const
221 {
222 return myValue;
223 }
224
225 Standard_Integer Index() const
226 {
227 return myMathIndex;
228 }
229
230 bool operator>(const SolInfo& theOther) const
231 {
232 return myValue > theOther.myValue;
233 }
234
235 bool operator<(const SolInfo& theOther) const
236 {
237 return myValue < theOther.myValue;
238 }
239
240 bool operator==(const SolInfo& theOther) const
241 {
242 return myValue == theOther.myValue;
243 }
244
245 Standard_Real& ChangeValue()
246 {
247 return myValue;
248 }
249
250private:
251 Standard_Integer myMathIndex;
252 Standard_Real myValue;
253};
254
74d80fb9 255static
256void BoundedArc (const TheArc& A,
257 const Handle(TheTopolTool)& Domain,
258 const Standard_Real Pdeb,
259 const Standard_Real Pfin,
260 TheFunction& Func,
261 IntStart_SequenceOfPathPoint& pnt,
262 IntStart_SequenceOfSegment& seg,
263 const Standard_Real TolBoundary,
264 const Standard_Real TolTangency,
265 Standard_Boolean& Arcsol,
266 const Standard_Boolean RecheckOnRegularity)
267{
f542b7bb 268 // Recherche des points solutions et des bouts d arc solution sur un arc donne.
269 // On utilise la fonction math_FunctionAllRoots. Ne convient donc que pour
270 // des arcs ayant un point debut et un point de fin (intervalle ferme de
271 // parametrage).
74d80fb9 272
9369e98a 273 Standard_Integer i, Nbi = 0, Nbp = 0;
74d80fb9 274
275 gp_Pnt ptdeb,ptfin;
295cb053 276 Standard_Real pardeb = 0., parfin = 0.;
74d80fb9 277 Standard_Integer ideb,ifin,range,ranged,rangef;
f542b7bb 278
74d80fb9 279 // Creer l echantillonage (math_FunctionSample ou classe heritant)
280 // Appel a math_FunctionAllRoots
281
74d80fb9 282 //@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@
283 //@@@ La Tolerance est asociee a l arc ( Incoherence avec le cheminement )
284 //@@@ ( EpsX ~ 1e-5 et ResolutionU et V ~ 1e-9 )
285 //@@@ le vertex trouve ici n'est pas retrouve comme point d arret d une
286 //@@@ ligne de cheminement
287 //@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@
3306fdd9 288 Standard_Real EpsX = 1.e-10;
74d80fb9 289 //@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@
290 //@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@
291 //@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@
f542b7bb 292
293 // Standard_Integer NbEchant = TheSOBTool::NbSamplesOnArc(A);
74d80fb9 294 Standard_Integer NbEchant = Func.NbSamples();
9369e98a 295 if(NbEchant<100) NbEchant = 100; //-- lbr le 22 Avril 96
296 //-- Toujours des pbs
297
74d80fb9 298 //-- Modif 24 Aout 93 -----------------------------
299 Standard_Real nTolTangency = TolTangency;
300 if((Pfin - Pdeb) < (TolTangency*10.0)) {
301 nTolTangency=(Pfin-Pdeb)*0.1;
302 }
303 if(EpsX>(nTolTangency+nTolTangency)) {
304 EpsX = nTolTangency * 0.1;
305 }
9369e98a 306
74d80fb9 307 //--------------------------------------------------
308 //-- Plante avec un edge avec 2 Samples
309 //-- dont les extremites son solutions (f=0)
310 //-- et ou la derivee est nulle
311 //-- Exemple : un segment diametre d une sphere
312 //-- if(NbEchant<3) NbEchant = 3; //-- lbr le 19 Avril 95
313 //--------------------------------------------------
314 Standard_Real para=0,dist,maxdist;
f542b7bb 315
74d80fb9 316 //-------------------------------------------------------------- REJECTIONS le 15 oct 98
317 Standard_Boolean Rejection=Standard_True;
318 Standard_Real maxdr,maxr,minr,ur,dur;
319 minr=RealLast();
320 maxr=-minr;
321 maxdr=-minr;
322 dur=(Pfin-Pdeb)*0.2;
323 for(i=1,ur=Pdeb;i<=6;i++) {
324 Standard_Real F,D;
325 if(Func.Values(ur,F,D)) {
326 Standard_Real lminr,lmaxr;
327 if(D<0.0) D=-D;
328 D*=dur+dur;
329 if(D>maxdr) maxdr=D;
330 lminr=F-D;
331 lmaxr=F+D;
332 if(lminr<minr) minr=lminr;
333 if(lmaxr>maxr) maxr=lmaxr;
334 if(minr<0.0 && maxr>0.0) {
335 Rejection=Standard_False;
78fdb3d3 336 break;
74d80fb9 337 }
338 }
339 ur+=dur;
340 }
78fdb3d3 341 if(Rejection)
342 {
343 dur=0.001+maxdr+(maxr-minr)*0.1;
344 minr-=dur;
345 maxr+=dur;
346 if(minr<0.0 && maxr>0.0) {
347 Rejection=Standard_False;
348 }
74d80fb9 349 }
350
351 Arcsol=Standard_False;
352
9369e98a 353 if(Rejection==Standard_False)
354 {
355 const IntSurf_Quadric& aQuadric = Func.Quadric();
356 GeomAbs_SurfaceType TypeQuad = aQuadric.TypeQuadric();
357
358 IntCurveSurface_HInter IntCS;
359 Standard_Boolean IsIntCSdone = Standard_False;
360 TColStd_SequenceOfReal Params;
361
362#if (defined(_MSC_VER) && (_MSC_VER < 1600))
363 std::auto_ptr<math_FunctionAllRoots> pSol;
364#else
365 std::unique_ptr<math_FunctionAllRoots> pSol;
366#endif
367
74d80fb9 368 math_FunctionSample Echant(Pdeb,Pfin,NbEchant);
f542b7bb 369
74d80fb9 370 Standard_Boolean aelargir=Standard_True;
371 //modified by NIZNHY-PKV Thu Apr 12 09:25:19 2001 f
372 //
373 //maxdist = 100.0*TolBoundary;
374 maxdist = TolBoundary+TolTangency;
375 //
376 //modified by NIZNHY-PKV Thu Apr 12 09:25:23 2001 t
377 for(i=1; i<=NbEchant && aelargir;i++) {
378 Standard_Real u = Echant.GetParameter(i);
379 if(Func.Value(u,dist)) {
380 if(dist>maxdist || -dist>maxdist) {
381 aelargir=Standard_False;
382 }
383 }
384 }
385 if(!(aelargir && maxdist<0.01)) {
386 maxdist = TolBoundary;
387 }
f542b7bb 388
9369e98a 389 if (TypeQuad != GeomAbs_OtherSurface) //intersection of boundary curve and quadric surface
390 {
391 //Exact solution
c22b52d6 392 Handle(Adaptor3d_Surface) aSurf = Func.Surface();
9369e98a 393 Adaptor3d_CurveOnSurface ConS(A, aSurf);
394 GeomAbs_CurveType TypeConS = ConS.GetType();
395#ifdef OCCT_DEBUG
396 Handle(Geom_Curve) CurveConS;
397 switch(TypeConS)
398 {
399 case GeomAbs_Line:
400 {
401 CurveConS = new Geom_Line(ConS.Line());
402 break;
403 }
404 case GeomAbs_Circle:
405 {
406 CurveConS = new Geom_Circle(ConS.Circle());
407 break;
408 }
409 case GeomAbs_Ellipse:
410 {
411 CurveConS = new Geom_Ellipse(ConS.Ellipse());
412 break;
413 }
414 case GeomAbs_Hyperbola:
415 {
416 CurveConS = new Geom_Hyperbola(ConS.Hyperbola());
417 break;
418 }
419 case GeomAbs_Parabola:
420 {
421 CurveConS = new Geom_Parabola(ConS.Parabola());
422 break;
423 }
424 case GeomAbs_BezierCurve:
425 {
426 CurveConS = ConS.Bezier();
427 break;
428 }
429 case GeomAbs_BSplineCurve:
430 {
431 CurveConS = ConS.BSpline();
432 break;
433 }
434 default:
435 {
436 Standard_Real MaxDeviation, AverageDeviation;
437 GeomLib::BuildCurve3d(1.e-5, ConS, ConS.FirstParameter(), ConS.LastParameter(),
438 CurveConS, MaxDeviation, AverageDeviation);
439 break;
440 }
441 }
442#endif
c22b52d6 443 Handle(Adaptor3d_CurveOnSurface) HConS = new Adaptor3d_CurveOnSurface(ConS);
9369e98a 444 Handle(Geom_Surface) QuadSurf;
445 switch (TypeQuad)
446 {
447 case GeomAbs_Plane:
448 {
449 QuadSurf = new Geom_Plane(aQuadric.Plane());
450 break;
451 }
452 case GeomAbs_Cylinder:
453 {
454 QuadSurf = new Geom_CylindricalSurface(aQuadric.Cylinder());
455 break;
456 }
457 case GeomAbs_Cone:
458 {
459 QuadSurf = new Geom_ConicalSurface(aQuadric.Cone());
460 break;
461 }
462 case GeomAbs_Sphere:
463 {
464 QuadSurf = new Geom_SphericalSurface(aQuadric.Sphere());
465 break;
466 }
467 case GeomAbs_Torus:
468 {
469 QuadSurf = new Geom_ToroidalSurface(aQuadric.Torus());
470 break;
471 }
472 default:
473 break;
474 }
c22b52d6 475 Handle(GeomAdaptor_Surface) GAHsurf = new GeomAdaptor_Surface(QuadSurf);
9369e98a 476
477 if ((TypeConS == GeomAbs_Line ||
478 TypeConS == GeomAbs_Circle ||
479 TypeConS == GeomAbs_Ellipse ||
480 TypeConS == GeomAbs_Parabola ||
481 TypeConS == GeomAbs_Hyperbola) &&
482 TypeQuad != GeomAbs_Torus &&
483 !IsDegenerated(HConS) &&
484 !IsDegenerated(aQuadric))
485 {
486 //exact intersection for only canonic curves and real quadric surfaces
487 IntCS.Perform(HConS, GAHsurf);
488 }
489
490 IsIntCSdone = IntCS.IsDone();
491 if (IsIntCSdone)
492 {
493 Nbp = IntCS.NbPoints();
494 Nbi = IntCS.NbSegments();
495 }
496 //If we have not got intersection, it may be touch with some tolerance,
497 //need to be checked
498 if (Nbp == 0 && Nbi == 0)
499 IsIntCSdone = Standard_False;
500
501 } //if (TypeQuad != GeomAbs_OtherSurface) - intersection of boundary curve and quadric surface
502
503 if (!IsIntCSdone)
504 {
505 pSol.reset(new math_FunctionAllRoots(Func,Echant,EpsX,maxdist,maxdist)); //-- TolBoundary,nTolTangency);
506
507 if (!pSol->IsDone()) {throw Standard_Failure();}
508
509 Nbp=pSol->NbPoints();
510 }
db005e48 511 //
74d80fb9 512 //jgv: build solution on the whole boundary
513 if (RecheckOnRegularity && Nbp > 0 && IsRegularity(A, Domain))
514 {
515 //Standard_Real theTol = Domain->MaxTolerance(A);
516 //theTol += theTol;
517 Standard_Real theTol = 5.e-4;
518 math_FunctionAllRoots SolAgain(Func,Echant,EpsX,theTol,theTol); //-- TolBoundary,nTolTangency);
f542b7bb 519
9775fa61 520 if (!SolAgain.IsDone()) {throw Standard_Failure();}
f542b7bb 521
74d80fb9 522 Standard_Integer Nbi_again = SolAgain.NbIntervals();
f542b7bb 523
74d80fb9 524 if (Nbi_again > 0)
525 {
526 Standard_Integer NbSamples = 10;
527 Standard_Real delta = (Pfin - Pdeb)/NbSamples;
528 Standard_Real GlobalTol = theTol*10;
529 Standard_Boolean SolOnBoundary = Standard_True;
530 for (i = 0; i <= NbSamples; i++)
531 {
532 Standard_Real aParam = Pdeb + i*delta;
533 Standard_Real aValue;
534 Func.Value(aParam, aValue);
535 if (Abs(aValue) > GlobalTol)
536 {
537 SolOnBoundary = Standard_False;
538 break;
539 }
540 }
541
542 if (SolOnBoundary)
543 {
544 for (i = 1; i <= Nbi_again; i++)
545 {
546 IntStart_TheSegment newseg;
547 newseg.SetValue(A);
548 // Recuperer point debut et fin, et leur parametre.
549 SolAgain.GetInterval(i,pardeb,parfin);
f542b7bb 550
74d80fb9 551 if (Abs(pardeb - Pdeb) <= Precision::PConfusion())
552 pardeb = Pdeb;
553 if (Abs(parfin - Pfin) <= Precision::PConfusion())
554 parfin = Pfin;
f542b7bb 555
74d80fb9 556 SolAgain.GetIntervalState(i,ideb,ifin);
f542b7bb 557
74d80fb9 558 //-- cout<<" Debug : IntStart_SearchOnBoundaries_1.gxx : i= "<<i<<" ParDeb:"<<pardeb<<" ParFin:"<<parfin<<endl;
f542b7bb 559
74d80fb9 560 ptdeb=Func.Valpoint(ideb);
561 ptfin=Func.Valpoint(ifin);
f542b7bb 562
74d80fb9 563 PointProcess(ptdeb,pardeb,A,Domain,pnt,theTol,ranged);
564 newseg.SetLimitPoint(pnt.Value(ranged),Standard_True);
565 PointProcess(ptfin,parfin,A,Domain,pnt,theTol,rangef);
566 newseg.SetLimitPoint(pnt.Value(rangef),Standard_False);
567 seg.Append(newseg);
568 }
569 Arcsol=Standard_True;
570 return;
571 }
572 }
9369e98a 573 } //if (RecheckOnRegularity && Nbp > 0 && IsRegularity(A, Domain))
74d80fb9 574 ////////////////////////////////////////////
f542b7bb 575
74d80fb9 576 //-- detection du cas ou la fonction est quasi tangente et que les
577 //-- zeros sont quasi confondus.
578 //-- Dans ce cas on prend le point "milieu"
579 //-- On suppose que les solutions sont triees.
580
74d80fb9 581 if(Nbp) {
f542b7bb 582 NCollection_Array1<SolInfo> aSI(1, Nbp);
583
584 for(i=1;i<=Nbp;i++)
585 {
9369e98a 586 if (IsIntCSdone)
587 aSI(i).Init(IntCS, i);
588 else
589 aSI(i).Init(*pSol, i);
74d80fb9 590 }
f542b7bb 591
592 std::sort(aSI.begin(), aSI.end());
593
74d80fb9 594 //modified by NIZNHY-PKV Wed Mar 21 18:34:18 2001 f
595 //////////////////////////////////////////////////////////
596 // The treatment of the situation when line(arc) that is
597 // tangent to cylinder(domain).
598 // We should have only one solution i.e Nbp=1. Ok?
599 // But we have 2,3,.. solutions. That is wrong ersult.
600 // The TreatLC(...) function is dedicated to solve the pb.
601 // PKV Fri Mar 23 12:17:29 2001
f542b7bb 602
9369e98a 603 Standard_Integer ip = TreatLC (A, Domain, aQuadric, TolBoundary, pnt);
74d80fb9 604 if (ip) {
f542b7bb 605 //////////////////////////////////////////////////////////
606 //modified by NIZNHY-PKV Wed Mar 21 18:34:23 2001 t
607 //
608 // Using of old usual way proposed by Laurent
609 //
610 for(i=1;i<Nbp;i++) {
611 Standard_Real parap1 = aSI(i + 1).Value();
612 para = aSI(i).Value();
613
614 Standard_Real param=(para+parap1)*0.5;
615 Standard_Real ym;
616 if(Func.Value(param,ym)) {
617 if(Abs(ym)<maxdist) {
618 // Modified by skv - Tue Aug 31 12:13:51 2004 OCC569 Begin
619 // Consider this interval as tangent one. Treat it to find
620 // parameter with the lowest function value.
621
622 // Compute the number of nodes.
623 Standard_Real aTol = TolBoundary*1000.0;
624 if(aTol > 0.001)
625 aTol = 0.001;
626
627 // fix floating point exception 569, chl-922-e9
628 parap1 = (Abs(parap1) < 1.e9) ? parap1 : ((parap1 >= 0.) ? 1.e9 : -1.e9);
629 para = (Abs(para) < 1.e9) ? para : ((para >= 0.) ? 1.e9 : -1.e9);
630
631 Standard_Integer aNbNodes = RealToInt(Ceiling((parap1 - para)/aTol));
632
633 Standard_Real aVal = RealLast();
634 //Standard_Integer aNbNodes = 23;
635 Standard_Real aDelta = (parap1 - para)/(aNbNodes + 1.);
636 Standard_Integer ii;
637 Standard_Real aCurPar;
638 Standard_Real aCurVal;
639
640 for (ii = 0; ii <= aNbNodes + 1; ii++) {
641 aCurPar = (ii < aNbNodes + 1) ? para + ii*aDelta : parap1;
642
643 if (Func.Value(aCurPar, aCurVal)) {
644 //if (aCurVal < aVal) {
645 if (Abs(aCurVal) < aVal) {
646 //aVal = aCurVal;
647 aVal = Abs(aCurVal);
648 param = aCurPar;
649 }
650 }
651 }
652 // Modified by skv - Tue Aug 31 12:13:51 2004 OCC569 End
653 aSI(i).ChangeValue() = Pdeb - 1;
654 aSI(i + 1).ChangeValue() = param;
74d80fb9 655 }
74d80fb9 656 }
657 }
f542b7bb 658
659 for (i=1; i<=Nbp; i++) {
660 para = aSI(i).Value();
661 if((para-Pdeb)<EpsX || (Pfin-para)<EpsX)
662 continue;
663
664 if(!Func.Value(para,dist))
665 continue;
666
667 dist = Abs(dist);
668
669 Standard_Integer anIndx = -1;
9369e98a 670 //const Standard_Real aParam = Sol->GetPoint(aSI(i).Index());
671 const Standard_Real aParam = aSI(i).Value();
f542b7bb 672 if (dist < maxdist)
673 {
9369e98a 674 if (!IsIntCSdone &&
675 (Abs(aParam - Pdeb) <= Precision::PConfusion() || Abs(aParam - Pfin) <= Precision::PConfusion()))
74d80fb9 676 {
9369e98a 677 anIndx = pSol->GetPointState(aSI(i).Index());
74d80fb9 678 }
f542b7bb 679 }
680
681 gp_Pnt aPnt(anIndx < 0 ? Func.LastComputedPoint() : Func.Valpoint(anIndx));
74d80fb9 682
f542b7bb 683 if (dist > 0.1*Precision::Confusion())
684 {
685 //Precise found points. It results in following:
686 // 1. Make the vertex nearer to the intersection line
687 // (see description to issue #27252 in order to
688 // understand necessity).
689 // 2. Merge two near vertices to single point.
690
691 //All members in TabSol array has already been sorted in increase order.
692 //Now, we limit precise boundaries in order to avoid changing this order.
693 const Standard_Real aFPar = (i == 1) ? Pdeb : (para + aSI(i - 1).Value()) / 2.0;
694 const Standard_Real aLPar = (i == Nbp) ? Pfin : (para + aSI(i + 1).Value()) / 2.0;
695
696 MinFunction aNewFunc(Func);
697 math_BrentMinimum aMin(Precision::Confusion());
698
699 aMin.Perform(aNewFunc, aFPar, para, aLPar);
700 if(aMin.IsDone())
701 {
702 para = aMin.Location();
703 const gp_Pnt2d aP2d(A->Value(para));
704 aPnt = Func.Surface()->Value(aP2d.X(), aP2d.Y());
705 }
74d80fb9 706 }
f542b7bb 707
708 PointProcess(aPnt, para, A, Domain, pnt, TolBoundary, range);
74d80fb9 709 }
f542b7bb 710 }// end of if(ip)
74d80fb9 711 } // end of if(Nbp)
712
713 // Pour chaque intervalle trouve faire
714 // Traiter les extremites comme des points
715 // Ajouter intervalle dans la liste des segments
f542b7bb 716
9369e98a 717 if (!IsIntCSdone)
718 Nbi = pSol->NbIntervals();
74d80fb9 719
720 if (!RecheckOnRegularity && Nbp) {
721 //--cout<<" Debug : IntStart_SearchOnBoundaries_1.gxx :Nbp>0 0 <- Nbi "<<Nbi<<endl;
722 Nbi=0;
723 }
724
725 //-- cout<<" Debug : IntStart_SearchOnBoundaries_1.gxx : Nbi : "<<Nbi<<endl;
f542b7bb 726
74d80fb9 727 for (i=1; i<=Nbi; i++) {
728 IntStart_TheSegment newseg;
729 newseg.SetValue(A);
730 // Recuperer point debut et fin, et leur parametre.
9369e98a 731 if (IsIntCSdone)
732 {
733 IntCurveSurface_IntersectionSegment IntSeg = IntCS.Segment(i);
734 IntCurveSurface_IntersectionPoint End1 = IntSeg.FirstPoint();
735 IntCurveSurface_IntersectionPoint End2 = IntSeg.SecondPoint();
736 pardeb = End1.W();
737 parfin = End2.W();
738 ptdeb = End1.Pnt();
739 ptfin = End2.Pnt();
740 }
741 else
742 {
743 pSol->GetInterval(i,pardeb,parfin);
744 pSol->GetIntervalState(i,ideb,ifin);
74d80fb9 745
9369e98a 746 //-- cout<<" Debug : IntStart_SearchOnBoundaries_1.gxx : i= "<<i<<" ParDeb:"<<pardeb<<" ParFin:"<<parfin<<endl;
747
748 ptdeb=Func.Valpoint(ideb);
749 ptfin=Func.Valpoint(ifin);
750 }
f542b7bb 751
74d80fb9 752 PointProcess(ptdeb,pardeb,A,Domain,pnt,TolBoundary,ranged);
753 newseg.SetLimitPoint(pnt.Value(ranged),Standard_True);
754 PointProcess(ptfin,parfin,A,Domain,pnt,TolBoundary,rangef);
755 newseg.SetLimitPoint(pnt.Value(rangef),Standard_False);
756 seg.Append(newseg);
757 }
758
759 if (Nbi==1) {
f542b7bb 760 if((Abs(pardeb - Pdeb) < Precision::PConfusion()) &&
761 (Abs(parfin - Pfin) < Precision::PConfusion()))
191478a5 762 {
74d80fb9 763 Arcsol=Standard_True;
764 }
765 }
766 }
767}
768
769//=======================================================================
770//function : ComputeBoundsfromInfinite
771//purpose :
772//=======================================================================
773// - PROVISIONAL - TEMPORARY - NOT GOOD - NYI - TO DO
774// - Temporary - temporary - not good - nyi - to do
775void ComputeBoundsfromInfinite(TheFunction& Func,
776 Standard_Real& PDeb,
777 Standard_Real& PFin,
778 Standard_Integer& NbEchant)
779{
780
781 // - We are looking for parameters for start and end of the arc (2d curve)
782 // - Infinity, a way to intersect the quadric with a portion of arc
783 // - Finished.
784 //
785 // - The quadric is a plane, a cylinder, a cone and a sphere.
786 // - Idea: We take any point on the arc and the fact grow
d94fa32e 787 // - Terminals to the signed distance function values or is likely
74d80fb9 788 // - S cancel.
789 //
790 // - WARNING: The following calculations provide a very estimated coarse parameters.
791 // - This avoids the raises and allows a case of Boxes
792 // - Inifinies walk. It will take this code
793 // - With curve surface intersections.
794
7c4e9501 795 NbEchant = 100;
74d80fb9 796
797 Standard_Real U0 = 0.0;
798 Standard_Real dU = 0.001;
799 Standard_Real Dist0,Dist1;
800
801 Func.Value(U0 , Dist0);
802 Func.Value(U0+dU, Dist1);
803 Standard_Real dDist = Dist1 - Dist0;
804 if(dDist) {
805 U0 -= dU*Dist0 / dDist;
806 PDeb = PFin = U0;
807 Standard_Real Umin = U0 - 1e5;
808 Func.Value(Umin , Dist0);
809 Func.Value(Umin+dU, Dist1);
810 dDist = Dist1-Dist0;
811 if(dDist) {
812 Umin -= dU*Dist0 / dDist;
813 }
814 else {
815 Umin-=10.0;
816 }
817 Standard_Real Umax = U0 + 1e8;
818 Func.Value(Umax , Dist0);
819 Func.Value(Umax+dU, Dist1);
820 dDist = Dist1-Dist0;
821 if(dDist) {
822 Umax -= dU*Dist0 / dDist;
823 }
824 else {
825 Umax+=10.0;
826 }
827 if(Umin>U0) { Umin=U0-10.0; }
828 if(Umax<U0) { Umax=U0+10.0; }
829
7c4e9501 830 PFin = Umax + 10. * (Umax - Umin);
831 PDeb = Umin - 10. * (Umax - Umin);
74d80fb9 832 }
833 else {
834 //-- Possibilite de Arc totalement inclu ds Quad
835 PDeb = 1e10;
836 PFin = -1e10;
837 }
838}
839
74d80fb9 840//=======================================================================
841//function : PointProcess
842//purpose :
843//=======================================================================
844void PointProcess (const gp_Pnt& Pt,
845 const Standard_Real Para,
846 const TheArc& A,
847 const Handle(TheTopolTool)& Domain,
848 IntStart_SequenceOfPathPoint& pnt,
849 const Standard_Real Tol,
850 Standard_Integer& Range)
851{
852
853// Check to see if a solution point is coincident with a vertex.
854// If confused, you should find this vertex in the list of
855// Start. It then returns the position of this point in the list pnt.
856// Otherwise, add the point in the list.
857
858 Standard_Integer k;
859 Standard_Boolean found,goon;
860 Standard_Real dist,toler;
861
862 Standard_Integer Nbsol = pnt.Length();
863 TheVertex vtx;
864 IntStart_ThePathPoint ptsol;
865
866 Domain->Initialize(A);
867 Domain->InitVertexIterator();
868 found = Standard_False;
869 goon = Domain->MoreVertex();
870 while (goon) {
871 vtx = Domain->Vertex();
872 dist= Abs(Para-TheSOBTool::Parameter(vtx,A));
873 toler = TheSOBTool::Tolerance(vtx,A);
0797d9d3 874#ifdef OCCT_DEBUG
74d80fb9 875 if(toler>0.1) {
04232180 876 std::cout<<"IntStart_SearchOnBoundaries_1.gxx : ** WARNING ** Tol Vertex="<<toler<<std::endl;
877 std::cout<<" Ou Edge degenere Ou Kro pointu"<<std::endl;
74d80fb9 878 if(toler>10000) toler=1e-7;
879 }
880#endif
881
882 if (dist <= toler) {
883 // Locate the vertex in the list of solutions
884 k=1;
885 found = (k>Nbsol);
886 while (!found) {
887 ptsol = pnt.Value(k);
888 if (!ptsol.IsNew()) {
889 //jag 940608 if (ptsol.Vertex() == vtx && ptsol.Arc() == A) {
890 if (Domain->Identical(ptsol.Vertex(),vtx) &&
891 ptsol.Arc() == A &&
892 Abs(ptsol.Parameter()-Para) <= toler) {
893 found=Standard_True;
894 }
895 else {
896 k=k+1;
897 found=(k>Nbsol);
898 }
899 }
900 else {
901 k=k+1;
902 found=(k>Nbsol);
903 }
904 }
905 if (k<=Nbsol) { // We find the vertex
906 Range = k;
907 }
908 else { // Otherwise
909 ptsol.SetValue(Pt,Tol,vtx,A,Para);
910 pnt.Append(ptsol);
911 Range = pnt.Length();
912 }
913 found = Standard_True;
914 goon = Standard_False;
915 }
916 else {
917 Domain->NextVertex();
918 goon = Domain->MoreVertex();
919 }
920 }
921
922 if (!found) { // No one is falling on a vertex
923 //jgv: do not add segment's extremities if they already exist
924 Standard_Boolean found_internal = Standard_False;
925 for (k = 1; k <= pnt.Length(); k++)
926 {
927 ptsol = pnt.Value(k);
928 if (ptsol.Arc() != A ||
929 !ptsol.IsNew()) //vertex
930 continue;
931 if (Abs(ptsol.Parameter()-Para) <= Precision::PConfusion())
932 {
933 found_internal = Standard_True;
934 Range = k;
935 }
936 }
937 /////////////////////////////////////////////////////////////
938
939 if (!found_internal)
940 {
941 Standard_Real TOL=Tol;
942 TOL*=1000.0;
db005e48 943 //if(TOL>0.001) TOL=0.001;
944 if(TOL>0.005) TOL=0.005; //#24643
74d80fb9 945
946 ptsol.SetValue(Pt,TOL,A,Para);
947 pnt.Append(ptsol);
948 Range = pnt.Length();
949 }
950 }
951}
952
953//=======================================================================
954//function : IsRegularity
955//purpose :
956//=======================================================================
957Standard_Boolean IsRegularity(const TheArc& /*A*/,
958 const Handle(TheTopolTool)& aDomain)
959{
960 Standard_Address anEAddress=aDomain->Edge();
961 if (anEAddress==NULL) {
962 return Standard_False;
963 }
964
965 TopoDS_Edge* anE=(TopoDS_Edge*)anEAddress;
966
967 return (BRep_Tool::HasContinuity(*anE));
968}
969
970//=======================================================================
971//function : TreatLC
972//purpose :
973//=======================================================================
974Standard_Integer TreatLC (const TheArc& A,
975 const Handle(TheTopolTool)& aDomain,
976 const IntSurf_Quadric& aQuadric,
977 const Standard_Real TolBoundary,
978 IntStart_SequenceOfPathPoint& pnt)
979{
980 Standard_Integer anExitCode=1, aNbExt;
981
982 Standard_Address anEAddress=aDomain->Edge();
983 if (anEAddress==NULL) {
984 return anExitCode;
985 }
986
987 TopoDS_Edge* anE=(TopoDS_Edge*)anEAddress;
988
989 if (BRep_Tool::Degenerated(*anE)) {
990 return anExitCode;
991 }
992
993 GeomAbs_CurveType aTypeE;
994 BRepAdaptor_Curve aBAC(*anE);
995 aTypeE=aBAC.GetType();
996
997 if (aTypeE!=GeomAbs_Line) {
998 return anExitCode;
999 }
1000
1001 GeomAbs_SurfaceType aTypeS;
1002 aTypeS=aQuadric.TypeQuadric();
1003
1004 if (aTypeS!=GeomAbs_Cylinder) {
1005 return anExitCode;
1006 }
1007
96a95605 1008 Standard_Real f, l, U1f, U1l, U2f, U2l, UEgde, TOL, aDist, aR, aRRel, Tol;
74d80fb9 1009 Handle(Geom_Curve) aCEdge=BRep_Tool::Curve(*anE, f, l);
1010
1011 gp_Cylinder aCyl=aQuadric.Cylinder();
1012 const gp_Ax1& anAx1=aCyl.Axis();
1013 gp_Lin aLin(anAx1);
1014 Handle(Geom_Line) aCAxis=new Geom_Line (aLin);
1015 aR=aCyl.Radius();
1016
1017 U1f = aCAxis->FirstParameter();
1018 U1l = aCAxis->LastParameter();
1019
1020 U2f = aCEdge->FirstParameter();
1021 U2l = aCEdge->LastParameter();
1022
1023
1024 GeomAdaptor_Curve C1, C2;
1025
1026 C1.Load(aCAxis);
1027 C2.Load(aCEdge);
1028
1029 Tol = Precision::PConfusion();
1030
1031 Extrema_ExtCC anExtCC(C1, C2, U1f, U1l, U2f, U2l, Tol, Tol);
1032
1033 aNbExt=anExtCC.NbExt();
1034 if (aNbExt!=1) {
1035 return anExitCode;
1036 }
1037
1038 gp_Pnt P1,PEdge;
1039 Extrema_POnCurv PC1, PC2;
1040
1041 anExtCC.Points(1, PC1, PC2);
1042
1043 P1 =PC1.Value();
1044 PEdge=PC2.Value();
1045
74d80fb9 1046 UEgde=PC2.Parameter();
1047
1048 aDist=PEdge.Distance(P1);
1049 aRRel=fabs(aDist-aR)/aR;
1050 if (aRRel > TolBoundary) {
1051 return anExitCode;
1052 }
1053
1054 if (UEgde < (f+TolBoundary) || UEgde > (l-TolBoundary)) {
1055 return anExitCode;
1056 }
1057 //
1058 // Do not wonder !
1059 // It was done as into PointProcess(...) function
1060 //printf("TreatLC()=> tangent line is found\n");
1061 TOL=1000.*TolBoundary;
1062 if(TOL>0.001) TOL=0.001;
1063
1064 IntStart_ThePathPoint ptsol;
1065 ptsol.SetValue(PEdge, TOL, A, UEgde);
1066 pnt.Append(ptsol);
1067
1068 anExitCode=0;
1069 return anExitCode;
1070
1071}
1072
1073
1074//=======================================================================
1075//function : IntStart_SearchOnBoundaries::IntStart_SearchOnBoundaries
1076//purpose :
1077//=======================================================================
1078IntStart_SearchOnBoundaries::IntStart_SearchOnBoundaries ()
d533dafb 1079: done(Standard_False),
1080 all(Standard_False)
74d80fb9 1081{
1082}
1083
1084//=======================================================================
1085//function : Perform
1086//purpose :
1087//=======================================================================
1088 void IntStart_SearchOnBoundaries::Perform (TheFunction& Func,
1089 const Handle(TheTopolTool)& Domain,
1090 const Standard_Real TolBoundary,
1091 const Standard_Real TolTangency,
1092 const Standard_Boolean RecheckOnRegularity)
1093{
1094
1095 done = Standard_False;
1096 spnt.Clear();
1097 sseg.Clear();
1098
1099 Standard_Boolean Arcsol;
1100 Standard_Real PDeb,PFin, prm, tol;
1101 Standard_Integer i, nbknown, nbfound,index;
1102 gp_Pnt pt;
1103
1104 Domain->Init();
1105
1106 if (Domain->More()) {
1107 all = Standard_True;
1108 }
1109 else {
1110 all = Standard_False;
1111 }
1112
1113 while (Domain->More()) {
1114 TheArc A = Domain->Value();
1115 if (!TheSOBTool::HasBeenSeen(A)) {
1116 Func.Set(A);
1117 FindVertex(A,Domain,Func,spnt,TolBoundary);
1118 TheSOBTool::Bounds(A,PDeb,PFin);
1119 if(Precision::IsNegativeInfinite(PDeb) ||
1120 Precision::IsPositiveInfinite(PFin)) {
7c4e9501 1121 Standard_Integer NbEchant;
1122 ComputeBoundsfromInfinite(Func,PDeb,PFin,NbEchant);
74d80fb9 1123 }
7c4e9501 1124 BoundedArc(A,Domain,PDeb,PFin,Func,spnt,sseg,TolBoundary,TolTangency,Arcsol,RecheckOnRegularity);
74d80fb9 1125 all = (all && Arcsol);
1126 }
1127
1128 else {
1129 // as it seems we'll never be here, because
1130 // TheSOBTool::HasBeenSeen(A) always returns FALSE
1131 nbfound = spnt.Length();
1132
1133 // On recupere les points connus
1134 nbknown = TheSOBTool::NbPoints(A);
1135 for (i=1; i<=nbknown; i++) {
1136 TheSOBTool::Value(A,i,pt,tol,prm);
1137 if (TheSOBTool::IsVertex(A,i)) {
1138 TheVertex vtx;
1139 TheSOBTool::Vertex(A,i,vtx);
1140 spnt.Append(IntStart_ThePathPoint(pt,tol,vtx,A,prm));
1141 }
1142 else {
1143 spnt.Append(IntStart_ThePathPoint(pt,tol,A,prm));
1144 }
1145 }
1146 // On recupere les arcs solutions
1147 nbknown = TheSOBTool::NbSegments(A);
1148 for (i=1; i<=nbknown; i++) {
1149 IntStart_TheSegment newseg;
1150 newseg.SetValue(A);
1151 if (TheSOBTool::HasFirstPoint(A,i,index)) {
1152 newseg.SetLimitPoint(spnt.Value(nbfound+index),Standard_True);
1153 }
1154 if (TheSOBTool::HasLastPoint(A,i,index)) {
1155 newseg.SetLimitPoint(spnt.Value(nbfound+index),Standard_False);
1156 }
1157 sseg.Append(newseg);
1158 }
1159 all = (all& TheSOBTool::IsAllSolution(A));
1160 }
1161 Domain->Next();
1162 }
1163 done = Standard_True;
1164}
Open CASCADE Technology repository