1 // Copyright (c) 1995-1999 Matra Datavision
2 // Copyright (c) 1999-2014 OPEN CASCADE SAS
4 // This file is part of Open CASCADE Technology software library.
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
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.
12 // Alternatively, this file may be used under the terms of Open CASCADE
13 // commercial license or contractual agreement.
16 #include <Adaptor3d_Surface.hxx>
17 #include <Adaptor3d_HSurfaceTool.hxx>
18 #include <Extrema_GenLocateExtPS.hxx>
19 #include <Geom_Surface.hxx>
22 #include <gp_Pnt2d.hxx>
23 #include <IntImp_ComputeTangence.hxx>
24 #include <IntSurf_LineOn2S.hxx>
25 #include <IntSurf_PntOn2S.hxx>
26 #include <IntWalk_PWalking.hxx>
27 #include <IntWalk_StatusDeflection.hxx>
28 #include <math_FunctionSetRoot.hxx>
29 #include <Precision.hxx>
30 #include <Standard_Failure.hxx>
31 #include <Standard_OutOfRange.hxx>
32 #include <StdFail_NotDone.hxx>
33 #include <TColgp_Array1OfPnt.hxx>
34 #include <TColStd_Array1OfReal.hxx>
36 //==================================================================================
37 // function : ComputePasInit
38 // purpose : estimate of max step : To avoid abrupt changes during change of isos
39 //==================================================================================
40 void IntWalk_PWalking::ComputePasInit(const Standard_Real theDeltaU1,
41 const Standard_Real theDeltaV1,
42 const Standard_Real theDeltaU2,
43 const Standard_Real theDeltaV2)
45 const Standard_Real aRangePart = 0.01;
46 const Standard_Real Increment = 2.0*pasMax;
47 const Handle(Adaptor3d_Surface)&
48 Caro1 = myIntersectionOn2S.Function().AuxillarSurface1();
49 const Handle(Adaptor3d_Surface)&
50 Caro2 = myIntersectionOn2S.Function().AuxillarSurface2();
52 const Standard_Real aDeltaU1=Abs(UM1-Um1);
53 const Standard_Real aDeltaV1=Abs(VM1-Vm1);
54 const Standard_Real aDeltaU2=Abs(UM2-Um2);
55 const Standard_Real aDeltaV2=Abs(VM2-Vm2);
57 //-- limit the reduction of uv box estimate to 0.01 natural box
58 //-- theDeltaU1 : On box of Inter
59 //-- aDeltaU1 : On parametric space
60 if(!Precision::IsInfinite(aDeltaU1))
61 pasuv[0]=Max(Increment*Max(theDeltaU1, aRangePart*aDeltaU1), pasuv[0]);
63 pasuv[0]=Max(Increment*theDeltaU1, pasuv[0]);
65 if(!Precision::IsInfinite(aDeltaV1))
66 pasuv[1]=Max(Increment*Max(theDeltaV1, aRangePart*aDeltaV1), pasuv[1]);
68 pasuv[1]=Max(Increment*theDeltaV1, pasuv[1]);
70 if(!Precision::IsInfinite(aDeltaU2))
71 pasuv[2]=Max(Increment*Max(theDeltaU2, aRangePart*aDeltaU2), pasuv[2]);
73 pasuv[2]=Max(Increment*theDeltaU2, pasuv[2]);
75 if(!Precision::IsInfinite(aDeltaV2))
76 pasuv[3]=Max(Increment*Max(theDeltaV2, aRangePart*aDeltaV2), pasuv[3]);
78 pasuv[3]=Max(Increment*theDeltaV2, pasuv[3]);
80 const Standard_Real ResoU1tol = Adaptor3d_HSurfaceTool::UResolution(Caro1, tolconf);
81 const Standard_Real ResoV1tol = Adaptor3d_HSurfaceTool::VResolution(Caro1, tolconf);
82 const Standard_Real ResoU2tol = Adaptor3d_HSurfaceTool::UResolution(Caro2, tolconf);
83 const Standard_Real ResoV2tol = Adaptor3d_HSurfaceTool::VResolution(Caro2, tolconf);
85 myStepMin[0] = Max(myStepMin[0], 2.0*ResoU1tol);
86 myStepMin[1] = Max(myStepMin[1], 2.0*ResoV1tol);
87 myStepMin[2] = Max(myStepMin[2], 2.0*ResoU2tol);
88 myStepMin[3] = Max(myStepMin[3], 2.0*ResoV2tol);
90 for(Standard_Integer i = 0; i < 4; i++)
92 pasuv[i]=Max(myStepMin[i], pasuv[i]);
96 //=======================================================================
97 //function : IsParallel
98 //purpose : Checks if theLine is parallel of some boundary of given
99 // surface (it is determined by theCheckSurf1 flag).
100 // Parallelism assumes small oscillations (swing is less or
101 // equal than theToler).
102 // Small lines (if first and last parameters in the Surface
103 // are almost equal) are classified as parallel (as same as
104 // any point can be considered as parallel of any line).
105 //=======================================================================
106 static void IsParallel(const Handle(IntSurf_LineOn2S)& theLine,
107 const Standard_Boolean theCheckSurf1,
108 const Standard_Real theToler,
109 Standard_Boolean& theIsUparallel,
110 Standard_Boolean& theIsVparallel)
112 const Standard_Integer aNbPointsMAX = 23;
114 theIsUparallel = theIsVparallel = Standard_True;
116 const Standard_Integer aNbLinePnts = theLine->NbPoints();
117 Standard_Integer aNbPoints = aNbLinePnts;
118 if(aNbPoints > aNbPointsMAX)
120 aNbPoints = aNbPointsMAX;
122 else if(aNbPoints < 3)
124 //Here we cannot estimate parallelism.
125 //Do all same as for small lines
129 Standard_Real aStep = IntToReal(theLine->NbPoints()) / aNbPoints;
130 Standard_Real aNPoint = 1.0;
132 Standard_Real aUmin = RealLast(), aUmax = RealFirst(), aVmin = RealLast(), aVmax = RealFirst();
133 for(Standard_Integer aNum = 1; aNum <= aNbPoints; aNum++, aNPoint += aStep)
135 // Fix possible "out of parameter" case.
136 if (aNPoint > aNbLinePnts)
137 aNPoint = aNbLinePnts;
141 theLine->Value(RealToInt(aNPoint)).ParametersOnS1(u, v);
143 theLine->Value(RealToInt(aNPoint)).ParametersOnS2(u, v);
158 theIsVparallel = ((aUmax - aUmin) < theToler);
159 theIsUparallel = ((aVmax - aVmin) < theToler);
162 //=======================================================================
163 //function : AdjustToDomain
164 //purpose : Returns TRUE if theP has been changed (i.e. initial value
165 // was out of the domain)
166 //=======================================================================
167 static Standard_Boolean AdjustToDomain(const Standard_Integer theNbElem,
168 Standard_Real* theParam,
169 const Standard_Real* const theLowBorder,
170 const Standard_Real* const theUppBorder)
172 Standard_Boolean aRetVal = Standard_False;
173 for (Standard_Integer i = 0; i < theNbElem; i++)
175 if ((theParam[i] - theLowBorder[i]) < -Precision::PConfusion())
177 theParam[i] = theLowBorder[i];
178 aRetVal = Standard_True;
181 if ((theParam[i] - theUppBorder[i]) > Precision::PConfusion())
183 theParam[i] = theUppBorder[i];
184 aRetVal = Standard_True;
191 //==================================================================================
192 // function : IntWalk_PWalking::IntWalk_PWalking
194 //==================================================================================
195 IntWalk_PWalking::IntWalk_PWalking(const Handle(Adaptor3d_Surface)& Caro1,
196 const Handle(Adaptor3d_Surface)& Caro2,
197 const Standard_Real TolTangency,
198 const Standard_Real Epsilon,
199 const Standard_Real Deflection,
200 const Standard_Real Increment )
204 close(Standard_False),
205 tgfirst(Standard_False),
206 tglast(Standard_False),
211 myTolTang(TolTangency),
213 previoustg(Standard_False),
214 myIntersectionOn2S(Caro1,Caro2,TolTangency),
215 STATIC_BLOCAGE_SUR_PAS_TROP_GRAND(0),
216 STATIC_PRECEDENT_INFLEXION(0)
218 Standard_Real KELARG=20.;
220 pasMax=Increment*0.2; //-- June 25 99 after problems with precision
221 Um1 = Adaptor3d_HSurfaceTool::FirstUParameter(Caro1);
222 Vm1 = Adaptor3d_HSurfaceTool::FirstVParameter(Caro1);
223 UM1 = Adaptor3d_HSurfaceTool::LastUParameter(Caro1);
224 VM1 = Adaptor3d_HSurfaceTool::LastVParameter(Caro1);
226 Um2 = Adaptor3d_HSurfaceTool::FirstUParameter(Caro2);
227 Vm2 = Adaptor3d_HSurfaceTool::FirstVParameter(Caro2);
228 UM2 = Adaptor3d_HSurfaceTool::LastUParameter(Caro2);
229 VM2 = Adaptor3d_HSurfaceTool::LastVParameter(Caro2);
231 ResoU1 = Adaptor3d_HSurfaceTool::UResolution(Caro1,Precision::Confusion());
232 ResoV1 = Adaptor3d_HSurfaceTool::VResolution(Caro1,Precision::Confusion());
234 ResoU2 = Adaptor3d_HSurfaceTool::UResolution(Caro2,Precision::Confusion());
235 ResoV2 = Adaptor3d_HSurfaceTool::VResolution(Caro2,Precision::Confusion());
237 Standard_Real NEWRESO;
238 Standard_Real MAXVAL;
239 Standard_Real MAXVAL2;
241 MAXVAL = Abs(Um1); MAXVAL2 = Abs(UM1);
242 if(MAXVAL2 > MAXVAL) MAXVAL = MAXVAL2;
243 NEWRESO = ResoU1 * MAXVAL ;
244 if(NEWRESO > ResoU1 &&NEWRESO<10) { ResoU1 = NEWRESO; }
247 MAXVAL = Abs(Um2); MAXVAL2 = Abs(UM2);
248 if(MAXVAL2 > MAXVAL) MAXVAL = MAXVAL2;
249 NEWRESO = ResoU2 * MAXVAL ;
250 if(NEWRESO > ResoU2 && NEWRESO<10) { ResoU2 = NEWRESO; }
253 MAXVAL = Abs(Vm1); MAXVAL2 = Abs(VM1);
254 if(MAXVAL2 > MAXVAL) MAXVAL = MAXVAL2;
255 NEWRESO = ResoV1 * MAXVAL ;
256 if(NEWRESO > ResoV1 && NEWRESO<10) { ResoV1 = NEWRESO; }
259 MAXVAL = Abs(Vm2); MAXVAL2 = Abs(VM2);
260 if(MAXVAL2 > MAXVAL) MAXVAL = MAXVAL2;
261 NEWRESO = ResoV2 * MAXVAL ;
262 if(NEWRESO > ResoV2 && NEWRESO<10) { ResoV2 = NEWRESO; }
264 pasuv[0]=pasMax*Abs(UM1-Um1);
265 pasuv[1]=pasMax*Abs(VM1-Vm1);
266 pasuv[2]=pasMax*Abs(UM2-Um2);
267 pasuv[3]=pasMax*Abs(VM2-Vm2);
269 if(ResoU1>0.0001*pasuv[0]) ResoU1=0.00001*pasuv[0];
270 if(ResoV1>0.0001*pasuv[1]) ResoV1=0.00001*pasuv[1];
271 if(ResoU2>0.0001*pasuv[2]) ResoU2=0.00001*pasuv[2];
272 if(ResoV2>0.0001*pasuv[3]) ResoV2=0.00001*pasuv[3];
275 if(Adaptor3d_HSurfaceTool::IsUPeriodic(Caro1)==Standard_False) {
276 //UM1+=KELARG*pasuv[0]; Um1-=KELARG*pasuv[0];
279 Standard_Real t = UM1-Um1;
280 if(t<Adaptor3d_HSurfaceTool::UPeriod(Caro1)) {
281 t=0.5*(Adaptor3d_HSurfaceTool::UPeriod(Caro1)-t);
282 t=(t>KELARG*pasuv[0])? KELARG*pasuv[0] : t;
287 if(Adaptor3d_HSurfaceTool::IsVPeriodic(Caro1)==Standard_False) {
288 //VM1+=KELARG*pasuv[1]; Vm1-=KELARG*pasuv[1];
291 Standard_Real t = VM1-Vm1;
292 if(t<Adaptor3d_HSurfaceTool::VPeriod(Caro1)) {
293 t=0.5*(Adaptor3d_HSurfaceTool::VPeriod(Caro1)-t);
294 t=(t>KELARG*pasuv[1])? KELARG*pasuv[1] : t;
299 if(Adaptor3d_HSurfaceTool::IsUPeriodic(Caro2)==Standard_False) {
300 //UM2+=KELARG*pasuv[2]; Um2-=KELARG*pasuv[2];
303 Standard_Real t = UM2-Um2;
304 if(t<Adaptor3d_HSurfaceTool::UPeriod(Caro2)) {
305 t=0.5*(Adaptor3d_HSurfaceTool::UPeriod(Caro2)-t);
306 t=(t>KELARG*pasuv[2])? KELARG*pasuv[2] : t;
311 if(Adaptor3d_HSurfaceTool::IsVPeriodic(Caro2)==Standard_False) {
312 //VM2+=KELARG*pasuv[3]; Vm2-=KELARG*pasuv[3];
315 Standard_Real t = VM2-Vm2;
316 if(t<Adaptor3d_HSurfaceTool::VPeriod(Caro2)) {
317 t=0.5*(Adaptor3d_HSurfaceTool::VPeriod(Caro2)-t);
318 t=(t>KELARG*pasuv[3])? KELARG*pasuv[3] : t;
323 myStepMin[0] = 100.0*ResoU1;
324 myStepMin[1] = 100.0*ResoV1;
325 myStepMin[2] = 100.0*ResoU2;
326 myStepMin[3] = 100.0*ResoV2;
328 //-- ComputePasInit(pasuv,Um1,UM1,Vm1,VM1,Um2,UM2,Vm2,VM2,Caro1,Caro2);
330 for (Standard_Integer i = 0; i<=3;i++) {
333 pasInit[i] = pasSav[i] = pasuv[i];
338 //==================================================================================
339 // function : IntWalk_PWalking
341 //==================================================================================
342 IntWalk_PWalking::IntWalk_PWalking(const Handle(Adaptor3d_Surface)& Caro1,
343 const Handle(Adaptor3d_Surface)& Caro2,
344 const Standard_Real TolTangency,
345 const Standard_Real Epsilon,
346 const Standard_Real Deflection,
347 const Standard_Real Increment,
348 const Standard_Real U1,
349 const Standard_Real V1,
350 const Standard_Real U2,
351 const Standard_Real V2)
355 close(Standard_False),
358 myTolTang(TolTangency),
360 myIntersectionOn2S(Caro1,Caro2,TolTangency),
361 STATIC_BLOCAGE_SUR_PAS_TROP_GRAND(0),
362 STATIC_PRECEDENT_INFLEXION(0)
364 Standard_Real KELARG=20.;
366 pasMax=Increment*0.2; //-- June 25 99 after problems with precision
368 Um1 = Adaptor3d_HSurfaceTool::FirstUParameter(Caro1);
369 Vm1 = Adaptor3d_HSurfaceTool::FirstVParameter(Caro1);
370 UM1 = Adaptor3d_HSurfaceTool::LastUParameter(Caro1);
371 VM1 = Adaptor3d_HSurfaceTool::LastVParameter(Caro1);
373 Um2 = Adaptor3d_HSurfaceTool::FirstUParameter(Caro2);
374 Vm2 = Adaptor3d_HSurfaceTool::FirstVParameter(Caro2);
375 UM2 = Adaptor3d_HSurfaceTool::LastUParameter(Caro2);
376 VM2 = Adaptor3d_HSurfaceTool::LastVParameter(Caro2);
378 ResoU1 = Adaptor3d_HSurfaceTool::UResolution(Caro1,Precision::Confusion());
379 ResoV1 = Adaptor3d_HSurfaceTool::VResolution(Caro1,Precision::Confusion());
381 ResoU2 = Adaptor3d_HSurfaceTool::UResolution(Caro2,Precision::Confusion());
382 ResoV2 = Adaptor3d_HSurfaceTool::VResolution(Caro2,Precision::Confusion());
384 Standard_Real NEWRESO, MAXVAL, MAXVAL2;
388 if(MAXVAL2 > MAXVAL) {
391 NEWRESO = ResoU1 * MAXVAL ;
392 if(NEWRESO > ResoU1) {
398 if(MAXVAL2 > MAXVAL){
401 NEWRESO = ResoU2 * MAXVAL ;
402 if(NEWRESO > ResoU2) {
408 if(MAXVAL2 > MAXVAL) {
411 NEWRESO = ResoV1 * MAXVAL ;
412 if(NEWRESO > ResoV1) {
418 if(MAXVAL2 > MAXVAL){
421 NEWRESO = ResoV2 * MAXVAL ;
422 if(NEWRESO > ResoV2) {
426 pasuv[0]=pasMax*Abs(UM1-Um1);
427 pasuv[1]=pasMax*Abs(VM1-Vm1);
428 pasuv[2]=pasMax*Abs(UM2-Um2);
429 pasuv[3]=pasMax*Abs(VM2-Vm2);
431 if(Adaptor3d_HSurfaceTool::IsUPeriodic(Caro1)==Standard_False) {
432 UM1+=KELARG*pasuv[0];
433 Um1-=KELARG*pasuv[0];
436 Standard_Real t = UM1-Um1;
437 if(t<Adaptor3d_HSurfaceTool::UPeriod(Caro1)) {
438 t=0.5*(Adaptor3d_HSurfaceTool::UPeriod(Caro1)-t);
439 t=(t>KELARG*pasuv[0])? KELARG*pasuv[0] : t;
445 if(Adaptor3d_HSurfaceTool::IsVPeriodic(Caro1)==Standard_False) {
446 VM1+=KELARG*pasuv[1];
447 Vm1-=KELARG*pasuv[1];
450 Standard_Real t = VM1-Vm1;
451 if(t<Adaptor3d_HSurfaceTool::VPeriod(Caro1)) {
452 t=0.5*(Adaptor3d_HSurfaceTool::VPeriod(Caro1)-t);
453 t=(t>KELARG*pasuv[1])? KELARG*pasuv[1] : t;
458 if(Adaptor3d_HSurfaceTool::IsUPeriodic(Caro2)==Standard_False) {
459 UM2+=KELARG*pasuv[2];
460 Um2-=KELARG*pasuv[2];
463 Standard_Real t = UM2-Um2;
464 if(t<Adaptor3d_HSurfaceTool::UPeriod(Caro2)) {
465 t=0.5*(Adaptor3d_HSurfaceTool::UPeriod(Caro2)-t);
466 t=(t>KELARG*pasuv[2])? KELARG*pasuv[2] : t;
472 if(Adaptor3d_HSurfaceTool::IsVPeriodic(Caro2)==Standard_False) {
473 VM2+=KELARG*pasuv[3];
474 Vm2-=KELARG*pasuv[3];
477 Standard_Real t = VM2-Vm2;
478 if(t<Adaptor3d_HSurfaceTool::VPeriod(Caro2)) {
479 t=0.5*(Adaptor3d_HSurfaceTool::VPeriod(Caro2)-t);
480 t=(t>KELARG*pasuv[3])? KELARG*pasuv[3] : t;
485 //-- ComputePasInit(pasuv,Um1,UM1,Vm1,VM1,Um2,UM2,Vm2,VM2,Caro1,Caro2);
487 for (Standard_Integer i = 0; i<=3;i++) {
488 pasInit[i] = pasSav[i] = pasuv[i];
491 if(ResoU1>0.0001*pasuv[0]) ResoU1=0.00001*pasuv[0];
492 if(ResoV1>0.0001*pasuv[1]) ResoV1=0.00001*pasuv[1];
493 if(ResoU2>0.0001*pasuv[2]) ResoU2=0.00001*pasuv[2];
494 if(ResoV2>0.0001*pasuv[3]) ResoV2=0.00001*pasuv[3];
496 myStepMin[0] = 100.0*ResoU1;
497 myStepMin[1] = 100.0*ResoV1;
498 myStepMin[2] = 100.0*ResoU2;
499 myStepMin[3] = 100.0*ResoV2;
502 TColStd_Array1OfReal Par(1,4);
510 //==================================================================================
511 // function : PerformFirstPoint
513 //==================================================================================
514 Standard_Boolean IntWalk_PWalking::PerformFirstPoint (const TColStd_Array1OfReal& ParDep,
515 IntSurf_PntOn2S& FirstPoint)
518 close = Standard_False;
521 TColStd_Array1OfReal Param(1,4);
523 for (i=1; i<=4; ++i) {
524 Param(i) = ParDep(i);
526 //-- calculate the first solution point
527 math_FunctionSetRoot Rsnld(myIntersectionOn2S.Function());
529 myIntersectionOn2S.Perform(Param,Rsnld);
530 if (!myIntersectionOn2S.IsDone()) {
531 return Standard_False;
534 if (myIntersectionOn2S.IsEmpty()) {
535 return Standard_False;
538 FirstPoint = myIntersectionOn2S.Point();
539 return Standard_True;
541 //==================================================================================
542 // function : Perform
544 //==================================================================================
545 void IntWalk_PWalking::Perform(const TColStd_Array1OfReal& ParDep)
547 Perform(ParDep,Um1,Vm1,Um2,Vm2,UM1,VM1,UM2,VM2);
550 //=======================================================================
551 //function : SQDistPointSurface
552 //purpose : Returns square distance between thePnt and theSurf.
553 // (theU0, theV0) is initial point for extrema
554 //=======================================================================
555 static Standard_Real SQDistPointSurface(const gp_Pnt &thePnt,
556 const Adaptor3d_Surface& theSurf,
557 const Standard_Real theU0,
558 const Standard_Real theV0)
560 Extrema_GenLocateExtPS aExtPS(theSurf);
561 aExtPS.Perform(thePnt, theU0, theV0);
566 return aExtPS.SquareDistance();
569 //==================================================================================
570 // function : IsTangentExtCheck
571 // purpose : Additional check if the surfaces are tangent.
572 // Checks if any point in one surface lie in another surface
573 // (with given tolerance)
574 //==================================================================================
575 static Standard_Boolean IsTangentExtCheck(const Handle(Adaptor3d_Surface)& theSurf1,
576 const Handle(Adaptor3d_Surface)& theSurf2,
577 const Standard_Real theU10,
578 const Standard_Real theV10,
579 const Standard_Real theU20,
580 const Standard_Real theV20,
581 const Standard_Real theToler,
582 const Standard_Real theArrStep[])
586 gp_Vec aDu1, aDv1, aDu2, aDv2;
587 theSurf1->D1(theU10, theV10, aPt, aDu1, aDv1);
588 theSurf2->D1(theU20, theV20, aPt, aDu2, aDv2);
590 const gp_Vec aN1(aDu1.Crossed(aDv1)),
591 aN2(aDu2.Crossed(aDv2));
592 const Standard_Real aDP = aN1.Dot(aN2),
593 aSQ1 = aN1.SquareMagnitude(),
594 aSQ2 = aN2.SquareMagnitude();
596 if((aSQ1 < RealSmall()) || (aSQ2 < RealSmall()))
597 return Standard_True; //Tangent
599 if(aDP*aDP < 0.9998*aSQ1*aSQ2)
600 {//cos(ang N1<->N2) < 0.9999
601 return Standard_False; //Not tangent
605 //For two faces (2^2 = 4)
606 const Standard_Real aSQToler = 4.0*theToler*theToler;
607 const Standard_Integer aNbItems = 4;
608 const Standard_Real aParUS1[aNbItems] = { theU10 + theArrStep[0],
609 theU10 - theArrStep[0],
611 const Standard_Real aParVS1[aNbItems] = { theV10, theV10,
612 theV10 + theArrStep[1],
613 theV10 - theArrStep[1]};
614 const Standard_Real aParUS2[aNbItems] = { theU20 + theArrStep[2],
615 theU20 - theArrStep[2],
617 const Standard_Real aParVS2[aNbItems] = { theV20, theV20,
618 theV20 + theArrStep[3],
619 theV20 - theArrStep[3]};
621 for(Standard_Integer i = 0; i < aNbItems; i++)
623 gp_Pnt aP(theSurf1->Value(aParUS1[i], aParVS1[i]));
624 const Standard_Real aSqDist = SQDistPointSurface(aP, *theSurf2, theU20, theV20);
625 if(aSqDist > aSQToler)
626 return Standard_False;
629 for(Standard_Integer i = 0; i < aNbItems; i++)
631 gp_Pnt aP(theSurf2->Value(aParUS2[i], aParVS2[i]));
632 const Standard_Real aSqDist = SQDistPointSurface(aP, *theSurf1, theU10, theV10);
633 if(aSqDist > aSQToler)
634 return Standard_False;
637 return Standard_True;
640 //==================================================================================
641 // function : Perform
643 //==================================================================================
644 void IntWalk_PWalking::Perform(const TColStd_Array1OfReal& ParDep,
645 const Standard_Real u1min,
646 const Standard_Real v1min,
647 const Standard_Real u2min,
648 const Standard_Real v2min,
649 const Standard_Real u1max,
650 const Standard_Real v1max,
651 const Standard_Real u2max,
652 const Standard_Real v2max)
654 const Standard_Real aSQDistMax = 1.0e-14;
657 Standard_Integer NbPasOKConseq=0;
658 TColStd_Array1OfReal Param(1,4);
659 IntImp_ConstIsoparametric ChoixIso;
662 done = Standard_False;
665 const Handle(Adaptor3d_Surface)& Caro1 =myIntersectionOn2S.Function().AuxillarSurface1();
666 const Handle(Adaptor3d_Surface)& Caro2 =myIntersectionOn2S.Function().AuxillarSurface2();
668 const Standard_Real UFirst1 = Adaptor3d_HSurfaceTool::FirstUParameter(Caro1);
669 const Standard_Real VFirst1 = Adaptor3d_HSurfaceTool::FirstVParameter(Caro1);
670 const Standard_Real ULast1 = Adaptor3d_HSurfaceTool::LastUParameter (Caro1);
671 const Standard_Real VLast1 = Adaptor3d_HSurfaceTool::LastVParameter (Caro1);
673 const Standard_Real UFirst2 = Adaptor3d_HSurfaceTool::FirstUParameter(Caro2);
674 const Standard_Real VFirst2 = Adaptor3d_HSurfaceTool::FirstVParameter(Caro2);
675 const Standard_Real ULast2 = Adaptor3d_HSurfaceTool::LastUParameter (Caro2);
676 const Standard_Real VLast2 = Adaptor3d_HSurfaceTool::LastVParameter (Caro2);
678 ComputePasInit(u1max - u1min,v1max - v1min,u2max - u2min,v2max - v2min);
680 for (Standard_Integer i=0; i<4; ++i)
687 pasInit[i] = pasSav[i] = pasuv[i];
690 line = new IntSurf_LineOn2S ();
692 for (Standard_Integer i=1; i<=4; ++i)
696 //-- reproduce steps uv connected to surfaces Caro1 and Caro2
697 //-- pasuv[] and pasSav[] are modified during the marching
698 for(Standard_Integer i = 0; i < 4; ++i)
700 pasSav[i] = pasuv[i] = pasInit[i];
703 //-- calculate the first solution point
704 math_FunctionSetRoot Rsnld(myIntersectionOn2S.Function());
706 ChoixIso = myIntersectionOn2S.Perform(Param,Rsnld);
707 if (!myIntersectionOn2S.IsDone())
713 if (myIntersectionOn2S.IsEmpty())
718 if(myIntersectionOn2S.IsTangent())
723 Standard_Boolean Arrive, DejaReparti;
724 const Standard_Integer RejectIndexMAX = 250000;
725 Standard_Integer IncKey, RejectIndex;
728 DejaReparti = Standard_False;
732 previousPoint = myIntersectionOn2S.Point();
733 previoustg = Standard_False;
734 previousd = myIntersectionOn2S.Direction();
735 previousd1 = myIntersectionOn2S.DirectionOnS1();
736 previousd2 = myIntersectionOn2S.DirectionOnS2();
739 firstd1 = previousd1;
740 firstd2 = previousd2;
741 tgfirst = tglast = Standard_False;
742 choixIsoSav = ChoixIso;
743 //------------------------------------------------------------
744 //-- Test if the first point of marching corresponds
745 //-- to a point on borders.
746 //-- In this case, DejaReparti is initialized as True
748 pf = previousPoint.Value();
749 Standard_Boolean bTestFirstPoint = Standard_True;
751 previousPoint.Parameters(Param(1),Param(2),Param(3),Param(4));
753 if(IsTangentExtCheck(Caro1, Caro2, Param(1), Param(2), Param(3), Param(4), myTolTang, pasuv))
756 AddAPoint(previousPoint);
758 IntWalk_StatusDeflection aStatus = IntWalk_OK, aPrevStatus = IntWalk_OK;
759 Standard_Boolean NoTestDeflection = Standard_False;
760 Standard_Real SvParam[4], f;
761 Standard_Integer LevelOfEmptyInmyIntersectionOn2S=0;
762 Standard_Integer LevelOfPointConfondu = 0;
763 Standard_Integer LevelOfIterWithoutAppend = -1;
766 const Standard_Real aTol[4] = { Epsilon(UM1 - Um1),
771 Standard_Integer aPrevNbPoints = line->NbPoints();
773 Arrive = Standard_False;
776 aPrevStatus = aStatus;
778 LevelOfIterWithoutAppend++;
779 if(LevelOfIterWithoutAppend>20)
781 Arrive = Standard_True;
785 RepartirOuDiviser(DejaReparti,ChoixIso,Arrive);
786 LevelOfIterWithoutAppend = 0;
792 case IntImp_UIsoparametricOnCaro1: f = Abs(previousd1.X()); break;
793 case IntImp_VIsoparametricOnCaro1: f = Abs(previousd1.Y()); break;
794 case IntImp_UIsoparametricOnCaro2: f = Abs(previousd2.X()); break;
795 case IntImp_VIsoparametricOnCaro2: f = Abs(previousd2.Y()); break;
803 previousPoint.Parameters(Param(1),Param(2),Param(3),Param(4));
806 Standard_Real aIncKey, aEps, dP1, dP2, dP3, dP4;
808 dP1 = sensCheminement * pasuv[0] * previousd1.X() /f;
809 dP2 = sensCheminement * pasuv[1] * previousd1.Y() /f;
810 dP3 = sensCheminement * pasuv[2] * previousd2.X() /f;
811 dP4 = sensCheminement * pasuv[3] * previousd2.Y() /f;
813 aIncKey=5.*(Standard_Real)IncKey;
815 if(ChoixIso == IntImp_UIsoparametricOnCaro1 && Abs(dP1) < aEps)
820 if(ChoixIso == IntImp_VIsoparametricOnCaro1 && Abs(dP2) < aEps)
825 if(ChoixIso == IntImp_UIsoparametricOnCaro2 && Abs(dP3) < aEps)
830 if(ChoixIso == IntImp_VIsoparametricOnCaro2 && Abs(dP4) < aEps)
840 //==========================
846 Standard_Integer aTryNumber = 0;
847 Standard_Boolean isBadPoint = Standard_False;
848 IntImp_ConstIsoparametric aBestIso = ChoixIso;
851 isBadPoint = Standard_False;
853 ChoixIso= myIntersectionOn2S.Perform(Param, Rsnld, aBestIso);
855 if (myIntersectionOn2S.IsDone() && !myIntersectionOn2S.IsEmpty())
857 //If we go along any surface boundary then it is possible
858 //to find "outboundaried" point.
859 //Nevertheless, if this deflection is quite small, we will be
860 //able to adjust this point to the boundary.
862 Standard_Real aNewPnt[4], anAbsParamDist[4];
863 myIntersectionOn2S.Point().Parameters(aNewPnt[0], aNewPnt[1], aNewPnt[2], aNewPnt[3]);
864 const Standard_Real aParMin[4] = {Um1, Vm1, Um2, Vm2};
865 const Standard_Real aParMax[4] = {UM1, VM1, UM2, VM2};
867 for(Standard_Integer i = 0; i < 4; i++)
869 if(Abs(aNewPnt[i] - aParMin[i]) < aTol[i])
870 aNewPnt[i] = aParMin[i];
871 else if(Abs(aNewPnt[i] - aParMax[i]) < aTol[i])
872 aNewPnt[i] = aParMax[i];
875 if (aNewPnt[0] < Um1 || aNewPnt[0] > UM1 ||
876 aNewPnt[1] < Vm1 || aNewPnt[1] > VM1 ||
877 aNewPnt[2] < Um2 || aNewPnt[2] > UM2 ||
878 aNewPnt[3] < Vm2 || aNewPnt[3] > VM2)
880 // Out of borders, process it later.
884 myIntersectionOn2S.ChangePoint().SetValue(aNewPnt[0],
889 anAbsParamDist[0] = Abs(Param(1) - dP1 - aNewPnt[0]);
890 anAbsParamDist[1] = Abs(Param(2) - dP2 - aNewPnt[1]);
891 anAbsParamDist[2] = Abs(Param(3) - dP3 - aNewPnt[2]);
892 anAbsParamDist[3] = Abs(Param(4) - dP4 - aNewPnt[3]);
893 if (anAbsParamDist[0] < ResoU1 &&
894 anAbsParamDist[1] < ResoV1 &&
895 anAbsParamDist[2] < ResoU2 &&
896 anAbsParamDist[3] < ResoV2 &&
897 aStatus != IntWalk_PasTropGrand)
899 isBadPoint = Standard_True;
900 aBestIso = IntImp_ConstIsoparametric((aBestIso + 1) % 4);
903 } while (isBadPoint && ++aTryNumber <= 4);
905 if (!myIntersectionOn2S.IsDone())
907 //end of line, division
908 Arrive = Standard_False;
913 RepartirOuDiviser(DejaReparti, ChoixIso, Arrive);
917 //== Calculation of exact point from Param(.) is possible
918 if (myIntersectionOn2S.IsEmpty())
920 Standard_Real u1,v1,u2,v2;
921 previousPoint.Parameters(u1,v1,u2,v2);
923 Arrive = Standard_False;
924 if(u1<UFirst1 || u1>ULast1)
926 Arrive=Standard_True;
929 if(u2<UFirst2 || u2>ULast2)
931 Arrive=Standard_True;
934 if(v1<VFirst1 || v1>VLast1)
936 Arrive=Standard_True;
939 if(v2<VFirst2 || v2>VLast2)
941 Arrive=Standard_True;
944 RepartirOuDiviser(DejaReparti,ChoixIso,Arrive);
945 LevelOfEmptyInmyIntersectionOn2S++;
947 if(LevelOfEmptyInmyIntersectionOn2S>10)
957 //============================================================
958 //== A point has been found : T E S T D E F L E C T I O N
959 //============================================================
962 NoTestDeflection = Standard_False;
966 if(--LevelOfEmptyInmyIntersectionOn2S<=0)
968 LevelOfEmptyInmyIntersectionOn2S=0;
969 if(LevelOfIterWithoutAppend < 10)
971 aStatus = TestDeflection(ChoixIso, aStatus);
975 if (aStatus != IntWalk_StepTooSmall)
977 // Bug #0029682 (Boolean intersection with
978 // fuzzy-option hangs).
979 // If aStatus == IntWalk_StepTooSmall then
980 // the counter "LevelOfIterWithoutAppend" will
981 // be nulified in the future if the step is smaller
982 // (see "case IntWalk_StepTooSmall:" below).
983 // Here, we forbid this nulification and thereby provide out from
994 //============================================================
995 //== T r a i t e m e n t s u r S t a t u s ==
996 //============================================================
997 if(LevelOfPointConfondu > 5)
999 aStatus = IntWalk_ArretSurPoint;
1000 LevelOfPointConfondu = 0;
1003 if(aStatus==IntWalk_OK)
1006 if(NbPasOKConseq >= 5)
1009 Standard_Boolean pastroppetit;
1014 pastroppetit=Standard_True;
1016 if(pasuv[0]<pasInit[0])
1018 t = (pasInit[0]-pasuv[0])*0.25;
1019 if(t>0.1*pasInit[0])
1025 pastroppetit=Standard_False;
1028 if(pasuv[1]<pasInit[1])
1030 t = (pasInit[1]-pasuv[1])*0.25;
1031 if(t>0.1*pasInit[1]) {
1036 pastroppetit=Standard_False;
1039 if(pasuv[2]<pasInit[2])
1041 t = (pasInit[2]-pasuv[2])*0.25;
1042 if(t>0.1*pasInit[2])
1048 pastroppetit=Standard_False;
1051 if(pasuv[3]<pasInit[3])
1053 t = (pasInit[3]-pasuv[3])*0.25;
1054 if(t>0.1*pasInit[3]) {
1058 pastroppetit=Standard_False;
1072 pastroppetit=Standard_False;
1076 while(pastroppetit);
1078 }//aStatus==IntWalk_OK
1083 switch(aStatus)//007
1085 case IntWalk_ArretSurPointPrecedent:
1087 Arrive = Standard_False;
1088 RepartirOuDiviser(DejaReparti, ChoixIso, Arrive);
1091 case IntWalk_PasTropGrand:
1093 Param(1)=SvParam[0];
1094 Param(2)=SvParam[1];
1095 Param(3)=SvParam[2];
1096 Param(4)=SvParam[3];
1098 if(LevelOfIterWithoutAppend > 5)
1100 for (Standard_Integer i = 0; i < 4; i++)
1102 if (pasSav[i] > pasInit[i])
1105 const Standard_Real aDelta = (pasInit[i]-pasSav[i])*0.25;
1107 if(aDelta > Epsilon(pasInit[i]))
1109 pasInit[i] -= aDelta;
1110 if ((aPrevStatus != IntWalk_StepTooSmall) &&
1111 (line->NbPoints() != aPrevNbPoints))
1113 LevelOfIterWithoutAppend = 0;
1116 aPrevNbPoints = line->NbPoints();
1123 case IntWalk_PointConfondu:
1125 LevelOfPointConfondu++;
1127 if(LevelOfPointConfondu>5)
1129 Standard_Boolean pastroppetit;
1133 pastroppetit=Standard_True;
1135 for(Standard_Integer i = 0; i < 4; i++)
1137 if(pasuv[i]<pasInit[i])
1139 pasuv[i]+=(pasInit[i]-pasuv[i])*0.25;
1140 pastroppetit=Standard_False;
1156 pastroppetit=Standard_False;
1160 while(pastroppetit);
1165 case IntWalk_StepTooSmall:
1167 Standard_Boolean hasStepBeenIncreased = Standard_False;
1169 for(Standard_Integer i = 0; i < 4; i++)
1171 const Standard_Real aNewStep = Min(1.5*pasuv[i], pasInit[i]);
1172 if(aNewStep > pasuv[i])
1174 pasuv[i] = aNewStep;
1175 hasStepBeenIncreased = Standard_True;
1179 if(hasStepBeenIncreased)
1181 Param(1)=SvParam[0];
1182 Param(2)=SvParam[1];
1183 Param(3)=SvParam[2];
1184 Param(4)=SvParam[3];
1186 // In order to avoid cyclic changes
1187 // (PasTropGrand --> Decrease step -->
1188 // StepTooSmall --> Increase step --> PasTropGrand...)
1189 // nullify LevelOfIterWithoutAppend only if the condition
1191 if ((aPrevStatus != IntWalk_PasTropGrand) &&
1192 (line->NbPoints() != aPrevNbPoints))
1194 LevelOfIterWithoutAppend = 0;
1197 aPrevNbPoints = line->NbPoints();
1202 Standard_FALLTHROUGH
1204 case IntWalk_ArretSurPoint://006
1206 //=======================================================
1207 //== Stop Test t : Frame on Param(.) ==
1208 //=======================================================
1210 Arrive = TestArret(DejaReparti,Param,ChoixIso);
1211 // JMB 30th December 1999.
1212 // Some statement below should not be put in comment because they are useful.
1213 // See grid CTO 909 A1 which infinitely loops
1214 if(Arrive==Standard_False && aStatus==IntWalk_ArretSurPoint)
1216 Arrive=Standard_True;
1218 std::cout << "IntWalk_PWalking_1.gxx: Problems with intersection"<<std::endl;
1224 NbPasOKConseq = -10;
1229 //=====================================================
1230 //== Param(.) is in the limits ==
1231 //== and does not end a closed line ==
1232 //=====================================================
1233 //== Check on the current point of myInters
1234 Standard_Boolean pointisvalid = Standard_False;
1236 Standard_Real u1,v1,u2,v2;
1237 myIntersectionOn2S.Point().Parameters(u1,v1,u2,v2);
1240 if(u1 <= UM1 && u2 <= UM2 && v1 <= VM1 &&
1241 v2 <= VM2 && u1 >= Um1 && u2 >= Um2 &&
1242 v1 >= Vm1 && v2 >= Vm2)
1244 pointisvalid=Standard_True;
1251 previousPoint = myIntersectionOn2S.Point();
1252 previoustg = myIntersectionOn2S.IsTangent();
1256 previousd = myIntersectionOn2S.Direction();
1257 previousd1 = myIntersectionOn2S.DirectionOnS1();
1258 previousd2 = myIntersectionOn2S.DirectionOnS2();
1260 //=====================================================
1261 //== Check on the previous Point
1263 Standard_Real u1,v1,u2,v2;
1264 previousPoint.Parameters(u1,v1,u2,v2);
1265 if( u1 <= UM1 && u2 <= UM2 && v1 <= VM1 &&
1266 v2 <= VM2 && u1 >= Um1 && u2 >= Um2 &&
1267 v1 >= Vm1 && v2 >= Vm2)
1269 pl = previousPoint.Value();
1272 if(pf.SquareDistance(pl) < aSQDistMax)
1282 bTestFirstPoint = Standard_False;
1286 AddAPoint(previousPoint);
1289 if(RejectIndex >= RejectIndexMAX)
1291 Arrive = Standard_True;
1296 LevelOfIterWithoutAppend = 0;
1300 //====================================================
1302 if (aStatus == IntWalk_ArretSurPoint)
1304 RepartirOuDiviser(DejaReparti,ChoixIso,Arrive);
1308 if (line->NbPoints() == 2)
1310 pasSav[0] = pasuv[0];
1311 pasSav[1] = pasuv[1];
1312 pasSav[2] = pasuv[2];
1313 pasSav[3] = pasuv[3];
1315 if ((aPrevStatus == IntWalk_PasTropGrand) &&
1316 (LevelOfIterWithoutAppend > 0))
1318 pasInit[0] = pasuv[0];
1319 pasInit[1] = pasuv[1];
1320 pasInit[2] = pasuv[2];
1321 pasInit[3] = pasuv[3];
1330 //================= la ligne est fermee ===============
1331 AddAPoint(line->Value(1)); //ligne fermee
1332 LevelOfIterWithoutAppend=0;
1336 //====================================================
1337 //== Param was not in the limits (was reframed)
1338 //====================================================
1339 Standard_Boolean bPrevNotTangent = !previoustg || !myIntersectionOn2S.IsTangent();
1341 IntImp_ConstIsoparametric SauvChoixIso = ChoixIso;
1342 ChoixIso = myIntersectionOn2S.Perform(Param,Rsnld,ChoixIso);
1344 if(!myIntersectionOn2S.IsEmpty()) //002
1346 // mutially outpasses in the square or intersection in corner
1348 if(TestArret(Standard_True,Param,ChoixIso))
1350 NbPasOKConseq = -10;
1351 ChoixIso = myIntersectionOn2S.Perform(Param,Rsnld,ChoixIso);
1353 if(!myIntersectionOn2S.IsEmpty())
1355 previousPoint = myIntersectionOn2S.Point();
1356 previoustg = myIntersectionOn2S.IsTangent();
1360 previousd = myIntersectionOn2S.Direction();
1361 previousd1 = myIntersectionOn2S.DirectionOnS1();
1362 previousd2 = myIntersectionOn2S.DirectionOnS2();
1365 pl = previousPoint.Value();
1369 if(pf.SquareDistance(pl) < aSQDistMax)
1379 bTestFirstPoint = Standard_False;
1383 AddAPoint(previousPoint);
1386 if(RejectIndex >= RejectIndexMAX)
1388 Arrive = Standard_True;
1393 LevelOfIterWithoutAppend=0;
1394 RepartirOuDiviser(DejaReparti,ChoixIso,Arrive);
1398 //fail framing divides the step
1399 Arrive = Standard_False;
1400 RepartirOuDiviser(DejaReparti,ChoixIso,Arrive);
1401 NoTestDeflection = Standard_True;
1402 ChoixIso = SauvChoixIso;
1407 // save the last point
1408 // to revert to it if the current point is out of bounds
1410 IntSurf_PntOn2S previousPointSave = previousPoint;
1411 Standard_Boolean previoustgSave = previoustg;
1412 gp_Dir previousdSave = previousd;
1413 gp_Dir2d previousd1Save = previousd1;
1414 gp_Dir2d previousd2Save = previousd2;
1416 previousPoint = myIntersectionOn2S.Point();
1417 previoustg = myIntersectionOn2S.IsTangent();
1418 Arrive = Standard_False;
1422 previousd = myIntersectionOn2S.Direction();
1423 previousd1 = myIntersectionOn2S.DirectionOnS1();
1424 previousd2 = myIntersectionOn2S.DirectionOnS2();
1427 //========================================
1428 //== Check on PreviousPoint @@
1431 Standard_Real u1,v1,u2,v2;
1432 previousPoint.Parameters(u1,v1,u2,v2);
1434 //To save initial 2d points
1435 gp_Pnt2d ParamPntOnS1(Param(1), Param(2));
1436 gp_Pnt2d ParamPntOnS2(Param(3), Param(4));
1438 ///////////////////////////
1446 Standard_Boolean bFlag1, bFlag2;
1447 Standard_Real aTol2D=1.e-11;
1449 bFlag1=u1 >= Um1-aTol2D && v1 >= Vm1-aTol2D && u1 <= UM1+aTol2D && v1 <= VM1+aTol2D;
1450 bFlag2=u2 >= Um2-aTol2D && v2 >= Vm2-aTol2D && u2 <= UM2+aTol2D && v2 <= VM2+aTol2D;
1451 if (bFlag1 && bFlag2)
1453 if (line->NbPoints() > 1)
1455 IntSurf_PntOn2S prevprevPoint = line->Value(line->NbPoints()-1);
1456 Standard_Real ppU1, ppV1, ppU2, ppV2;
1457 prevprevPoint.Parameters(ppU1, ppV1, ppU2, ppV2);
1458 Standard_Real pU1, pV1, pU2, pV2;
1459 previousPointSave.Parameters(pU1, pV1, pU2, pV2);
1460 gp_Vec2d V1onS1(gp_Pnt2d(ppU1, ppV1), gp_Pnt2d(pU1, pV1));
1461 gp_Vec2d V2onS1(gp_Pnt2d(pU1, pV1), gp_Pnt2d(u1, v1));
1462 gp_Vec2d V1onS2(gp_Pnt2d(ppU2, ppV2), gp_Pnt2d(pU2, pV2));
1463 gp_Vec2d V2onS2(gp_Pnt2d(pU2, pV2), gp_Pnt2d(u2, v2));
1465 const Standard_Real aDot1 = V1onS1 * V2onS1;
1466 const Standard_Real aDot2 = V1onS2 * V2onS2;
1468 if ((aDot1 < 0.0) || (aDot2 < 0.0))
1470 Arrive = Standard_True;
1475 if(u1 <= UM1 && u2 <= UM2 && v1 <= VM1 &&
1476 v2 <= VM2 && u1 >= Um1 && u2 >= Um2 &&
1477 v1 >= Vm1 && v2 >= Vm2) {
1480 pl = previousPoint.Value();
1484 if(pf.SquareDistance(pl) < aSQDistMax)
1495 bTestFirstPoint = Standard_False;
1499 //To avoid walking around the same point
1500 //in the tangent zone near a border
1504 //There are three consecutive points:
1505 //previousPointSave -> ParamPnt -> curPnt.
1507 Standard_Real prevU1, prevV1, prevU2, prevV2;
1508 previousPointSave.Parameters(prevU1, prevV1, prevU2, prevV2);
1509 gp_Pnt2d prevPntOnS1(prevU1, prevV1), prevPntOnS2(prevU2, prevV2);
1510 gp_Pnt2d curPntOnS1(u1, v1), curPntOnS2(u2, v2);
1511 gp_Vec2d PrevToParamOnS1(prevPntOnS1, ParamPntOnS1);
1512 gp_Vec2d PrevToCurOnS1(prevPntOnS1, curPntOnS1);
1513 gp_Vec2d PrevToParamOnS2(prevPntOnS2, ParamPntOnS2);
1514 gp_Vec2d PrevToCurOnS2(prevPntOnS2, curPntOnS2);
1515 Standard_Real MaxAngle = 3*M_PI/4;
1516 Standard_Real anAngleS1 = 0.0, anAngleS2 = 0.0;
1517 const Standard_Real aSQMParS1 = PrevToParamOnS1.SquareMagnitude();
1518 const Standard_Real aSQMParS2 = PrevToParamOnS2.SquareMagnitude();
1519 const Standard_Real aSQMCurS1 = PrevToCurOnS1.SquareMagnitude();
1520 const Standard_Real aSQMCurS2 = PrevToCurOnS2.SquareMagnitude();
1522 if(aSQMCurS1 < gp::Resolution())
1524 //We came back to the one of previous point.
1525 //Therefore, we must break;
1529 else if(aSQMParS1 < gp::Resolution())
1531 //We are walking along tangent zone.
1532 //It should be continued.
1537 anAngleS1 = Abs(PrevToParamOnS1.Angle(PrevToCurOnS1));
1540 if(aSQMCurS2 < gp::Resolution())
1542 //We came back to the one of previous point.
1543 //Therefore, we must break;
1547 else if(aSQMParS2 < gp::Resolution())
1549 //We are walking along tangent zone.
1550 //It should be continued;
1555 anAngleS2 = Abs(PrevToParamOnS2.Angle(PrevToCurOnS2));
1558 if ((anAngleS1 > MaxAngle) && (anAngleS2 > MaxAngle))
1560 Arrive = Standard_True;
1565 //Check singularity.
1566 //I.e. check if we are walking along direction, which does not
1567 //result in coming to any point (i.e. derivative
1568 //3D-intersection curve along this direction is equal to 0).
1569 //A sphere with direction {dU=1, dV=0} from point
1570 //(U=0, V=M_PI/2) can be considered as example for
1571 //this case (we cannot find another 3D-point if we go thus).
1573 //Direction chosen along 1st and 2nd surface correspondingly
1574 const gp_Vec2d aDirS1(prevPntOnS1, curPntOnS1),
1575 aDirS2(prevPntOnS2, curPntOnS2);
1578 gp_Vec aDuS1, aDvS1, aDuS2, aDvS2;
1580 myIntersectionOn2S.Function().AuxillarSurface1()->
1581 D1(curPntOnS1.X(), curPntOnS1.Y(), aPtemp, aDuS1, aDvS1);
1582 myIntersectionOn2S.Function().AuxillarSurface2()->
1583 D1(curPntOnS2.X(), curPntOnS2.Y(), aPtemp, aDuS2, aDvS2);
1585 //Derivative WLine along (it is vector-function indeed)
1587 //(https://en.wikipedia.org/wiki/Directional_derivative#Variation_using_only_direction_of_vector).
1588 //F1 - on the 1st surface, F2 - on the 2nd surface.
1589 //x, y, z - coordinates of derivative vector.
1590 const Standard_Real aF1x = aDuS1.X()*aDirS1.X() +
1591 aDvS1.X()*aDirS1.Y();
1592 const Standard_Real aF1y = aDuS1.Y()*aDirS1.X() +
1593 aDvS1.Y()*aDirS1.Y();
1594 const Standard_Real aF1z = aDuS1.Z()*aDirS1.X() +
1595 aDvS1.Z()*aDirS1.Y();
1596 const Standard_Real aF2x = aDuS2.X()*aDirS2.X() +
1597 aDvS2.X()*aDirS2.Y();
1598 const Standard_Real aF2y = aDuS2.Y()*aDirS2.X() +
1599 aDvS2.Y()*aDirS2.Y();
1600 const Standard_Real aF2z = aDuS2.Z()*aDirS2.X() +
1601 aDvS2.Z()*aDirS2.Y();
1603 const Standard_Real aF1 = aF1x*aF1x + aF1y*aF1y + aF1z*aF1z;
1604 const Standard_Real aF2 = aF2x*aF2x + aF2y*aF2y + aF2z*aF2z;
1606 if((aF1 < gp::Resolution()) && (aF2 < gp::Resolution()))
1608 //All derivative are equal to 0. Therefore, there is
1609 //no point in going along direction chosen.
1610 Arrive = Standard_True;
1614 }//if (previoustg) cond.
1616 ////////////////////////////////////////
1617 AddAPoint(previousPoint);
1620 if(RejectIndex >= RejectIndexMAX)
1622 Arrive = Standard_True;
1628 LevelOfIterWithoutAppend=0;
1629 Arrive = Standard_True;
1633 // revert to the last correctly calculated point
1634 previousPoint = previousPointSave;
1635 previoustg = previoustgSave;
1636 previousd = previousdSave;
1637 previousd1 = previousd1Save;
1638 previousd2 = previousd2Save;
1643 Standard_Boolean wasExtended = Standard_False;
1645 if(Arrive && myIntersectionOn2S.IsTangent() && bPrevNotTangent)
1647 if(ExtendLineInCommonZone(SauvChoixIso, DejaReparti))
1649 wasExtended = Standard_True;
1650 Arrive = Standard_False;
1651 ChoixIso = SauvChoixIso;
1655 RepartirOuDiviser(DejaReparti,ChoixIso,Arrive);
1658 myIntersectionOn2S.IsDone() && !myIntersectionOn2S.IsEmpty() &&
1659 myIntersectionOn2S.IsTangent() && bPrevNotTangent &&
1662 if(ExtendLineInCommonZone(SauvChoixIso, DejaReparti))
1664 wasExtended = Standard_True;
1665 Arrive = Standard_False;
1666 ChoixIso = SauvChoixIso;
1669 }//else !TestArret() $
1670 }//$$ end successful framing on border (!myIntersectionOn2S.IsEmpty())
1673 //echec framing on border; division of step
1674 Arrive = Standard_False;
1675 NoTestDeflection = Standard_True;
1676 RepartirOuDiviser(DejaReparti,ChoixIso,Arrive);
1678 }//$$$ end framing on border (!close)
1679 }//004 fin TestArret return Arrive = True
1680 } // 006case IntWalk_ArretSurPoint: end Processing Status = OK or ArretSurPoint
1681 } //007 switch(aStatus)
1682 } //008 end processing point (TEST DEFLECTION)
1683 } //009 end processing line (else if myIntersectionOn2S.IsDone())
1684 } //010 end if first departure point allows marching while (!Arrive)
1686 done = Standard_True;
1688 // ===========================================================================================================
1689 // function: ExtendLineInCommonZone
1690 // purpose: Extends already computed line inside tangent zone in the direction given by theChoixIso.
1691 // Returns Standard_True if the line was extended through tangent zone and the last computed point
1692 // is outside the tangent zone (but it is not put into the line). Otherwise returns Standard_False.
1693 // ===========================================================================================================
1694 Standard_Boolean IntWalk_PWalking::ExtendLineInCommonZone(const IntImp_ConstIsoparametric theChoixIso,
1695 const Standard_Boolean theDirectionFlag)
1698 const Handle(Adaptor3d_Surface)& Caro1 = myIntersectionOn2S.Function().AuxillarSurface1();
1699 const Handle(Adaptor3d_Surface)& Caro2 = myIntersectionOn2S.Function().AuxillarSurface2();
1701 const Standard_Real UFirst1 = Adaptor3d_HSurfaceTool::FirstUParameter(Caro1);
1702 const Standard_Real VFirst1 = Adaptor3d_HSurfaceTool::FirstVParameter(Caro1);
1703 const Standard_Real ULast1 = Adaptor3d_HSurfaceTool::LastUParameter(Caro1);
1704 const Standard_Real VLast1 = Adaptor3d_HSurfaceTool::LastVParameter(Caro1);
1706 const Standard_Real UFirst2 = Adaptor3d_HSurfaceTool::FirstUParameter(Caro2);
1707 const Standard_Real VFirst2 = Adaptor3d_HSurfaceTool::FirstVParameter(Caro2);
1708 const Standard_Real ULast2 = Adaptor3d_HSurfaceTool::LastUParameter(Caro2);
1709 const Standard_Real VLast2 = Adaptor3d_HSurfaceTool::LastVParameter(Caro2);
1711 Standard_Boolean bOutOfTangentZone = Standard_False;
1712 Standard_Boolean bStop = !myIntersectionOn2S.IsTangent();
1713 Standard_Integer dIncKey = 1;
1714 TColStd_Array1OfReal Param(1,4);
1715 IntWalk_StatusDeflection aStatus = IntWalk_OK;
1716 Standard_Integer nbIterWithoutAppend = 0;
1717 Standard_Integer nbEqualPoints = 0;
1718 Standard_Integer parit = 0;
1719 Standard_Integer uvit = 0;
1720 IntSurf_SequenceOfPntOn2S aSeqOfNewPoint;
1722 previousPoint.Parameters(Param(1), Param(2), Param(3), Param(4));
1724 if (Param(1) - UFirst1 < ResoU1)
1726 return bOutOfTangentZone;
1728 else if (Param(2) - VFirst1 < ResoV1)
1730 return bOutOfTangentZone;
1732 else if (Param(3) - UFirst2 < ResoU2)
1734 return bOutOfTangentZone;
1736 else if (Param(4) - VFirst2 < ResoV2)
1738 return bOutOfTangentZone;
1741 if (Param(1) - ULast1 > -ResoU1)
1743 return bOutOfTangentZone;
1745 else if (Param(2) - VLast1 > -ResoV1)
1747 return bOutOfTangentZone;
1749 else if (Param(3) - ULast2 > -ResoU2)
1751 return bOutOfTangentZone;
1753 else if (Param(4) - VLast2 > -ResoV2)
1755 return bOutOfTangentZone;
1759 nbIterWithoutAppend++;
1761 if((nbIterWithoutAppend > 20) || (nbEqualPoints > 20)) {
1763 std::cout<<"Infinite loop detected. Stop iterations (IntWalk_PWalking_1.gxx)" << std::endl;
1765 bStop = Standard_True;
1768 Standard_Real f = 0.;
1770 switch (theChoixIso)
1772 case IntImp_UIsoparametricOnCaro1: f = Abs(previousd1.X()); break;
1773 case IntImp_VIsoparametricOnCaro1: f = Abs(previousd1.Y()); break;
1774 case IntImp_UIsoparametricOnCaro2: f = Abs(previousd2.X()); break;
1775 case IntImp_VIsoparametricOnCaro2: f = Abs(previousd2.Y()); break;
1780 previousPoint.Parameters(Param(1),Param(2),Param(3),Param(4));
1782 Standard_Real dP1 = sensCheminement * pasuv[0] * previousd1.X() /f;
1783 Standard_Real dP2 = sensCheminement * pasuv[1] * previousd1.Y() /f;
1784 Standard_Real dP3 = sensCheminement * pasuv[2] * previousd2.X() /f;
1785 Standard_Real dP4 = sensCheminement * pasuv[3] * previousd2.Y() /f;
1787 if(theChoixIso == IntImp_UIsoparametricOnCaro1 && Abs(dP1) < 1.e-7) dP1 *= (5. * (Standard_Real)dIncKey);
1788 if(theChoixIso == IntImp_VIsoparametricOnCaro1 && Abs(dP2) < 1.e-7) dP2 *= (5. * (Standard_Real)dIncKey);
1789 if(theChoixIso == IntImp_UIsoparametricOnCaro2 && Abs(dP3) < 1.e-7) dP3 *= (5. * (Standard_Real)dIncKey);
1790 if(theChoixIso == IntImp_VIsoparametricOnCaro2 && Abs(dP4) < 1.e-7) dP4 *= (5. * (Standard_Real)dIncKey);
1796 Standard_Real SvParam[4];
1797 IntImp_ConstIsoparametric ChoixIso = theChoixIso;
1799 for(parit = 0; parit < 4; parit++) {
1800 SvParam[parit] = Param(parit+1);
1802 math_FunctionSetRoot Rsnld(myIntersectionOn2S.Function());
1803 ChoixIso = myIntersectionOn2S.Perform(Param,Rsnld, theChoixIso);
1805 if (!myIntersectionOn2S.IsDone()) {
1806 return bOutOfTangentZone;
1809 if (myIntersectionOn2S.IsEmpty()) {
1810 return bOutOfTangentZone;
1813 aStatus = TestDeflection(ChoixIso, aStatus);
1815 if(aStatus == IntWalk_OK) {
1817 for(uvit = 0; uvit < 4; uvit++) {
1818 if(pasuv[uvit] < pasInit[uvit]) {
1819 pasuv[uvit] = pasInit[uvit];
1825 case IntWalk_ArretSurPointPrecedent:
1827 bStop = Standard_True;
1828 bOutOfTangentZone = !myIntersectionOn2S.IsTangent();
1831 case IntWalk_PasTropGrand:
1833 for(parit = 0; parit < 4; parit++) {
1834 Param(parit+1) = SvParam[parit];
1836 Standard_Boolean bDecrease = Standard_False;
1838 for(uvit = 0; uvit < 4; uvit++) {
1839 if(pasSav[uvit] < pasInit[uvit]) {
1840 pasInit[uvit] -= (pasInit[uvit] - pasSav[uvit]) * 0.1;
1841 bDecrease = Standard_True;
1845 if(bDecrease) nbIterWithoutAppend--;
1848 case IntWalk_PointConfondu:
1850 for(uvit = 0; uvit < 4; uvit++) {
1851 if(pasuv[uvit] < pasInit[uvit]) {
1852 pasuv[uvit] += (pasInit[uvit] - pasuv[uvit]) * 0.1;
1858 case IntWalk_ArretSurPoint:
1861 bStop = TestArret(theDirectionFlag, Param, ChoixIso);
1866 Standard_Real u11,v11,u12,v12;
1867 myIntersectionOn2S.Point().Parameters(u11,v11,u12,v12);
1868 Standard_Real u21,v21,u22,v22;
1869 previousPoint.Parameters(u21,v21,u22,v22);
1871 if(((fabs(u11-u21) < ResoU1) && (fabs(v11-v21) < ResoV1)) ||
1872 ((fabs(u12-u22) < ResoU2) && (fabs(v12-v22) < ResoV2))) {
1881 bStop = bStop || !myIntersectionOn2S.IsTangent();
1882 bOutOfTangentZone = !myIntersectionOn2S.IsTangent();
1885 Standard_Boolean pointisvalid = Standard_False;
1886 Standard_Real u1,v1,u2,v2;
1887 myIntersectionOn2S.Point().Parameters(u1,v1,u2,v2);
1889 if(u1 <= UM1 && u2 <= UM2 && v1 <= VM1 &&
1890 v2 <= VM2 && u1 >= Um1 && u2 >= Um2 &&
1891 v1 >= Vm1 && v2 >= Vm2)
1892 pointisvalid = Standard_True;
1895 previousPoint = myIntersectionOn2S.Point();
1896 previoustg = myIntersectionOn2S.IsTangent();
1899 previousd = myIntersectionOn2S.Direction();
1900 previousd1 = myIntersectionOn2S.DirectionOnS1();
1901 previousd2 = myIntersectionOn2S.DirectionOnS2();
1903 Standard_Boolean bAddPoint = Standard_True;
1905 if(line->NbPoints() >= 1) {
1906 gp_Pnt pf = line->Value(1).Value();
1907 gp_Pnt pl = previousPoint.Value();
1909 if(pf.Distance(pl) < Precision::Confusion()) {
1911 if(dIncKey == 5000) return bOutOfTangentZone;
1912 else bAddPoint = Standard_False;
1917 aSeqOfNewPoint.Append(previousPoint);
1918 nbIterWithoutAppend = 0;
1922 if (line->NbPoints() == 2) {
1923 for(uvit = 0; uvit < 4; uvit++) {
1924 pasSav[uvit] = pasuv[uvit];
1928 if ( !pointisvalid ) {
1929 // decrease step if out of bounds
1930 // otherwise the same calculations will be
1931 // repeated several times
1932 if ( ( u1 > UM1 ) || ( u1 < Um1 ) )
1935 if ( ( v1 > VM1 ) || ( v1 < Vm1 ) )
1938 if ( ( u2 > UM2 ) || ( u2 < Um2 ) )
1941 if ( ( v2 > VM2 ) || ( v2 < Vm2 ) )
1946 if(close && (line->NbPoints() >= 1)) {
1948 if(!bOutOfTangentZone) {
1949 aSeqOfNewPoint.Append(line->Value(1)); // line end
1951 nbIterWithoutAppend = 0;
1954 ChoixIso = myIntersectionOn2S.Perform(Param, Rsnld, theChoixIso);
1956 if(myIntersectionOn2S.IsEmpty()) {
1957 bStop = Standard_True;// !myIntersectionOn2S.IsTangent();
1958 bOutOfTangentZone = Standard_False; // !myIntersectionOn2S.IsTangent();
1961 Standard_Boolean bAddPoint = Standard_True;
1962 Standard_Boolean pointisvalid = Standard_False;
1964 previousPoint = myIntersectionOn2S.Point();
1965 Standard_Real u1,v1,u2,v2;
1966 previousPoint.Parameters(u1,v1,u2,v2);
1968 if(u1 <= UM1 && u2 <= UM2 && v1 <= VM1 &&
1969 v2 <= VM2 && u1 >= Um1 && u2 >= Um2 &&
1970 v1 >= Vm1 && v2 >= Vm2)
1971 pointisvalid = Standard_True;
1975 if(line->NbPoints() >= 1) {
1976 gp_Pnt pf = line->Value(1).Value();
1977 gp_Pnt pl = previousPoint.Value();
1979 if(pf.Distance(pl) < Precision::Confusion()) {
1981 if(dIncKey == 5000) return bOutOfTangentZone;
1982 else bAddPoint = Standard_False;
1986 if(bAddPoint && !bOutOfTangentZone) {
1987 aSeqOfNewPoint.Append(previousPoint);
1988 nbIterWithoutAppend = 0;
2003 Standard_Boolean bExtendLine = Standard_False;
2004 Standard_Real u1 = 0., v1 = 0., u2 = 0., v2 = 0.;
2006 Standard_Integer pit = 0;
2008 for(pit = 0; !bExtendLine && (pit < 2); pit++) {
2010 previousPoint.Parameters(u1,v1,u2,v2);
2012 if(aSeqOfNewPoint.Length() > 0)
2013 aSeqOfNewPoint.Value(aSeqOfNewPoint.Length()).Parameters(u1,v1,u2,v2);
2018 if(((u1 - Um1) < ResoU1) ||
2019 ((UM1 - u1) < ResoU1) ||
2020 ((u2 - Um2) < ResoU2) ||
2021 ((UM2 - u2) < ResoU2) ||
2022 ((v1 - Vm1) < ResoV1) ||
2023 ((VM1 - v1) < ResoV1) ||
2024 ((v2 - Vm2) < ResoV2) ||
2025 ((VM2 - v2) < ResoV2))
2026 bExtendLine = Standard_True;
2030 // if(aStatus == IntWalk_OK || aStatus == IntWalk_ArretSurPoint) {
2031 if(aStatus == IntWalk_OK) {
2032 bExtendLine = Standard_True;
2034 if(aSeqOfNewPoint.Length() > 1) {
2035 TColStd_Array1OfReal FirstParams(0, 3), LastParams(0, 3), Resolutions(0, 3);
2036 Resolutions(0) = ResoU1; Resolutions(1) = ResoV1; Resolutions(2) = ResoU2; Resolutions(3) = ResoV2;
2038 aSeqOfNewPoint(1).Parameters(FirstParams.ChangeValue(0), FirstParams.ChangeValue(1),
2039 FirstParams.ChangeValue(2), FirstParams.ChangeValue(3));
2040 aSeqOfNewPoint(aSeqOfNewPoint.Length()).Parameters(LastParams.ChangeValue(0),
2041 LastParams.ChangeValue(1),
2042 LastParams.ChangeValue(2),
2043 LastParams.ChangeValue(3));
2044 Standard_Integer indexofiso = 0;
2046 if(theChoixIso == IntImp_UIsoparametricOnCaro1) indexofiso = 0;
2047 if(theChoixIso == IntImp_VIsoparametricOnCaro1) indexofiso = 1;
2048 if(theChoixIso == IntImp_UIsoparametricOnCaro2) indexofiso = 2;
2049 if(theChoixIso == IntImp_VIsoparametricOnCaro2) indexofiso = 3;
2051 Standard_Integer afirstindex = (indexofiso < 2) ? 0 : 2;
2052 gp_Vec2d aTangentZoneDir(gp_Pnt2d(FirstParams.Value(afirstindex), FirstParams.Value(afirstindex + 1)),
2053 gp_Pnt2d(LastParams.Value(afirstindex), LastParams.Value(afirstindex + 1)));
2055 gp_Dir2d anIsoDir(0, 1);
2057 if((indexofiso == 1) || (indexofiso == 3))
2058 anIsoDir = gp_Dir2d(1, 0);
2060 if(aTangentZoneDir.SquareMagnitude() > gp::Resolution()) {
2061 Standard_Real piquota = M_PI*0.25;
2063 if(fabs(aTangentZoneDir.Angle(anIsoDir)) > piquota) {
2064 Standard_Integer ii = 1, nextii = 2;
2066 Standard_Real asqresol = gp::Resolution();
2067 asqresol *= asqresol;
2070 aSeqOfNewPoint(ii).Parameters(FirstParams.ChangeValue(0), FirstParams.ChangeValue(1),
2071 FirstParams.ChangeValue(2), FirstParams.ChangeValue(3));
2072 aSeqOfNewPoint(ii + 1).Parameters(LastParams.ChangeValue(0), LastParams.ChangeValue(1),
2073 LastParams.ChangeValue(2), LastParams.ChangeValue(3));
2074 d1 = gp_Vec2d(gp_Pnt2d(FirstParams.Value(afirstindex),
2075 FirstParams.Value(afirstindex + 1)),
2076 gp_Pnt2d(LastParams.Value(afirstindex),
2077 LastParams.Value(afirstindex + 1)));
2080 while((d1.SquareMagnitude() < asqresol) &&
2081 (ii < aSeqOfNewPoint.Length()));
2085 while(nextii < aSeqOfNewPoint.Length()) {
2087 gp_Vec2d nextd1(0, 0);
2088 Standard_Integer jj = nextii;
2091 aSeqOfNewPoint(jj).Parameters(FirstParams.ChangeValue(0), FirstParams.ChangeValue(1),
2092 FirstParams.ChangeValue(2), FirstParams.ChangeValue(3));
2093 aSeqOfNewPoint(jj + 1).Parameters(LastParams.ChangeValue(0), LastParams.ChangeValue(1),
2094 LastParams.ChangeValue(2), LastParams.ChangeValue(3));
2095 nextd1 = gp_Vec2d(gp_Pnt2d(FirstParams.Value(afirstindex),
2096 FirstParams.Value(afirstindex + 1)),
2097 gp_Pnt2d(LastParams.Value(afirstindex),
2098 LastParams.Value(afirstindex + 1)));
2102 while((nextd1.SquareMagnitude() < asqresol) &&
2103 (jj < aSeqOfNewPoint.Length()));
2106 if(fabs(d1.Angle(nextd1)) > piquota) {
2107 bExtendLine = Standard_False;
2113 // end if(fabs(aTangentZoneDir.Angle(anIsoDir)
2120 return Standard_False;
2122 Standard_Integer i = 0;
2124 for(i = 1; i <= aSeqOfNewPoint.Length(); i++) {
2125 AddAPoint(aSeqOfNewPoint.Value(i));
2128 return bOutOfTangentZone;
2131 //=======================================================================
2132 //function : DistanceMinimizeByGradient
2136 // theInit should be initialized before function calling.
2137 //=======================================================================
2138 Standard_Boolean IntWalk_PWalking::
2139 DistanceMinimizeByGradient( const Handle(Adaptor3d_Surface)& theASurf1,
2140 const Handle(Adaptor3d_Surface)& theASurf2,
2141 TColStd_Array1OfReal& theInit,
2142 const Standard_Real* theStep0)
2144 const Standard_Integer aNbIterMAX = 60;
2145 const Standard_Real aTol = 1.0e-14;
2146 const Standard_Real aTolNul = 1.0 / Precision::Infinite();
2148 // I.e. if theU1 = 0.0 then Epsilon(theU1) = DBL_MIN (~1.0e-308).
2149 // Work with this number is impossible: there is a dangerous to
2150 // obtain Floating-point-overflow. Therefore, we limit this value.
2151 const Standard_Real aMinAddValU1 = Max(Epsilon(theInit(1)), aTolNul);
2152 const Standard_Real aMinAddValV1 = Max(Epsilon(theInit(2)), aTolNul);
2153 const Standard_Real aMinAddValU2 = Max(Epsilon(theInit(3)), aTolNul);
2154 const Standard_Real aMinAddValV2 = Max(Epsilon(theInit(4)), aTolNul);
2156 Handle(Geom_Surface) aS1, aS2;
2158 if (theASurf1->GetType() != GeomAbs_BezierSurface &&
2159 theASurf1->GetType() != GeomAbs_BSplineSurface)
2160 return Standard_True;
2161 if (theASurf2->GetType() != GeomAbs_BezierSurface &&
2162 theASurf2->GetType() != GeomAbs_BSplineSurface)
2163 return Standard_True;
2165 Standard_Boolean aStatus = Standard_False;
2168 gp_Vec aD1u, aD1v, aD2U, aD2V;
2170 theASurf1->D1(theInit(1), theInit(2), aP1, aD1u, aD1v);
2171 theASurf2->D1(theInit(3), theInit(4), aP2, aD2U, aD2V);
2173 Standard_Real aSQDistPrev = aP1.SquareDistance(aP2);
2175 gp_Vec aP12(aP1, aP2);
2177 Standard_Real aGradFu(-aP12.Dot(aD1u));
2178 Standard_Real aGradFv(-aP12.Dot(aD1v));
2179 Standard_Real aGradFU( aP12.Dot(aD2U));
2180 Standard_Real aGradFV( aP12.Dot(aD2V));
2182 Standard_Real aStepU1 = 1.0e-6, aStepV1 = 1.0e-6,
2183 aStepU2 = 1.0e-6, aStepV2 = 1.0e-6;
2187 aStepU1 = theStep0[0];
2188 aStepV1 = theStep0[1];
2189 aStepU2 = theStep0[2];
2190 aStepV2 = theStep0[3];
2193 Standard_Boolean flRepeat = Standard_True;
2194 Standard_Integer aNbIter = aNbIterMAX;
2198 Standard_Real anAdd = aGradFu*aStepU1;
2199 const Standard_Real aPARu = theInit(1) - Sign(Max(Abs(anAdd), aMinAddValU1), anAdd);
2201 anAdd = aGradFv*aStepV1;
2202 const Standard_Real aPARv = theInit(2) - Sign(Max(Abs(anAdd), aMinAddValV1), anAdd);
2204 anAdd = aGradFU*aStepU2;
2205 const Standard_Real aParU = theInit(3) - Sign(Max(Abs(anAdd), aMinAddValU2), anAdd);
2207 anAdd = aGradFV*aStepV2;
2208 const Standard_Real aParV = theInit(4) - Sign(Max(Abs(anAdd), aMinAddValV2), anAdd);
2212 theASurf1->D1(aPARu, aPARv, aPt1, aD1u, aD1v);
2213 theASurf2->D1(aParU, aParV, aPt2, aD2U, aD2V);
2215 Standard_Real aSQDist = aPt1.SquareDistance(aPt2);
2217 if(aSQDist < aSQDistPrev)
2219 aSQDistPrev = aSQDist;
2225 aStatus = aSQDistPrev < aTol;
2235 flRepeat = Standard_False;
2239 theASurf1->D1(theInit(1), theInit(2), aPt1, aD1u, aD1v);
2240 theASurf2->D1(theInit(3), theInit(4), aPt2, aD2U, aD2V);
2242 gp_Vec aPt12(aPt1, aPt2);
2243 aGradFu = -aPt12.Dot(aD1u);
2244 aGradFv = -aPt12.Dot(aD1v);
2245 aGradFU = aPt12.Dot(aD2U);
2246 aGradFV = aPt12.Dot(aD2V);
2250 aStepU1 = theStep0[0];
2251 aStepV1 = theStep0[1];
2252 aStepU2 = theStep0[2];
2253 aStepV2 = theStep0[3];
2257 aStepU1 = aStepV1 = aStepU2 = aStepV2 = 1.0e-6;
2266 //=======================================================================
2267 //function : DistanceMinimizeByExtrema
2271 // theP0, theU0 and theV0 parameters should be initialized
2272 // before the function calling.
2273 //=======================================================================
2274 Standard_Boolean IntWalk_PWalking::
2275 DistanceMinimizeByExtrema(const Handle(Adaptor3d_Surface)& theASurf,
2276 const gp_Pnt& theP0,
2277 Standard_Real& theU0,
2278 Standard_Real& theV0,
2279 const Standard_Real* theStep0)
2281 const Standard_Real aTol = 1.0e-14;
2283 gp_Vec aD1Su, aD1Sv, aD2Su, aD2Sv, aD2SuvTemp;
2284 Standard_Real aSQDistPrev = RealLast();
2285 Standard_Real aU = theU0, aV = theV0;
2287 const Standard_Real aStep0[2] = { theStep0 ? theStep0[0] : 1.0,
2288 theStep0 ? theStep0[1] : 1.0 };
2290 Standard_Integer aNbIter = 10;
2293 theASurf->D2(aU, aV, aPS, aD1Su, aD1Sv, aD2Su, aD2Sv, aD2SuvTemp);
2295 gp_Vec aVec(theP0, aPS);
2297 Standard_Real aSQDist = aVec.SquareMagnitude();
2299 if(aSQDist >= aSQDistPrev)
2302 aSQDistPrev = aSQDist;
2307 if(aSQDistPrev < aTol)
2311 const Standard_Real aF1 = aD1Su.Dot(aVec), aF2 = aD1Sv.Dot(aVec);
2314 const Standard_Real aDf1u = aD2Su.Dot(aVec) + aD1Su.Dot(aD1Su),
2315 aDf1v = aD2Su.Dot(aD1Sv),
2317 aDf2v = aD2Sv.Dot(aVec) + aD1Sv.Dot(aD1Sv);
2319 const Standard_Real aDet = aDf1u*aDf2v - aDf1v*aDf2u;
2320 aU -= aStep0[0]*(aDf2v*aF1 - aDf1v*aF2) / aDet;
2321 aV += aStep0[1]*(aDf2u*aF1 - aDf1u*aF2) / aDet;
2325 return (aSQDistPrev < aTol);
2328 //=======================================================================
2329 //function : HandleSingleSingularPoint
2331 //=======================================================================
2332 Standard_Boolean IntWalk_PWalking::HandleSingleSingularPoint(const Handle(Adaptor3d_Surface)& theASurf1,
2333 const Handle(Adaptor3d_Surface)& theASurf2,
2334 const Standard_Real the3DTol,
2335 TColStd_Array1OfReal &thePnt)
2337 // u1, v1, u2, v2 order is used.
2338 Standard_Real aLowBorder[4] = {theASurf1->FirstUParameter(),
2339 theASurf1->FirstVParameter(),
2340 theASurf2->FirstUParameter(),
2341 theASurf2->FirstVParameter()};
2342 Standard_Real aUppBorder[4] = {theASurf1->LastUParameter(),
2343 theASurf1->LastVParameter(),
2344 theASurf2->LastUParameter(),
2345 theASurf2->LastVParameter()};
2346 IntImp_ConstIsoparametric aLockedDir[4] = {IntImp_UIsoparametricOnCaro1,
2347 IntImp_VIsoparametricOnCaro1,
2348 IntImp_UIsoparametricOnCaro2,
2349 IntImp_VIsoparametricOnCaro2};
2351 // Create new intersector with new tolerance.
2352 IntWalk_TheInt2S anInt(theASurf1, theASurf2, the3DTol);
2353 math_FunctionSetRoot aRsnld(anInt.Function());
2355 for (Standard_Integer i = 1; i <= 4; ++i)
2357 if ( Abs(thePnt(i) - aLowBorder[i - 1]) < Precision::PConfusion() ||
2358 Abs(thePnt(i) - aUppBorder[i - 1]) < Precision::PConfusion())
2361 anInt.Perform(thePnt,aRsnld, aLockedDir[i - 1]);
2363 if (!anInt.IsDone())
2366 if (anInt.IsEmpty())
2369 Standard_Real aPars[4];
2370 anInt.Point().Parameters(aPars[0], aPars[1], aPars[2], aPars[3]);
2371 Handle(Adaptor3d_Surface) aSurfs[2] = { theASurf1, theASurf2 };
2373 Standard_Real aTol2 = the3DTol * the3DTol;
2378 for (k = 0; k < 2; ++k)
2380 Standard_Integer iu, iv;
2383 aSurfs[k]->D1(aPars[iu], aPars[iv], aPInt, aDU, aDV);
2384 Standard_Real aMod = aDU.Magnitude();
2385 if (aMod > Precision::Confusion())
2387 Standard_Real aTolU = the3DTol / aMod;
2388 if (Abs(aLowBorder[iu] - aPars[iu]) < aTolU)
2390 aP = aSurfs[k]->Value(aLowBorder[iu], aPars[iv]);
2391 if (aPInt.SquareDistance(aP) < aTol2)
2392 aPars[iu] = aLowBorder[iu];
2394 else if (Abs(aUppBorder[iu] - aPars[iu]) < aTolU)
2396 aP = aSurfs[k]->Value(aUppBorder[iu], aPars[iv]);
2397 if (aPInt.SquareDistance(aP) < aTol2)
2398 aPars[iu] = aUppBorder[iu];
2401 aMod = aDV.Magnitude();
2402 if (aMod > Precision::Confusion())
2404 Standard_Real aTolV = the3DTol / aMod;
2405 if (Abs(aLowBorder[iv] - aPars[iv]) < aTolV)
2407 aP = aSurfs[k]->Value(aPars[iu], aLowBorder[iv]);
2408 if (aPInt.SquareDistance(aP) < aTol2)
2409 aPars[iv] = aLowBorder[iv];
2411 else if (Abs(aUppBorder[iv] - aPars[iv]) < aTolV)
2413 aP = aSurfs[k]->Value(aPars[iu], aUppBorder[iv]);
2414 if (aPInt.SquareDistance(aP) < aTol2)
2415 aPars[iv] = aUppBorder[iv];
2422 for (j = 1; j <= 4; ++j)
2424 thePnt(j) = aPars[j - 1];
2426 Standard_Boolean isInDomain = Standard_True;
2427 for (j = 1; isInDomain && (j <= 4); ++j)
2429 if ((thePnt(j) - aLowBorder[j - 1] + Precision::PConfusion())*
2430 (thePnt(j) - aUppBorder[j - 1] - Precision::PConfusion()) > 0.0)
2432 isInDomain = Standard_False;
2437 return Standard_True;
2441 return Standard_False;
2444 //=======================================================================
2445 //function : SeekPointOnBoundary
2447 //=======================================================================
2448 Standard_Boolean IntWalk_PWalking::
2449 SeekPointOnBoundary(const Handle(Adaptor3d_Surface)& theASurf1,
2450 const Handle(Adaptor3d_Surface)& theASurf2,
2451 const Standard_Real theU1,
2452 const Standard_Real theV1,
2453 const Standard_Real theU2,
2454 const Standard_Real theV2,
2455 const Standard_Boolean isTheFirst)
2457 Standard_Boolean isOK = Standard_False;
2459 // u1, v1, u2, v2 order is used.
2460 const Standard_Real aLowBorder[4] = {theASurf1->FirstUParameter(),
2461 theASurf1->FirstVParameter(),
2462 theASurf2->FirstUParameter(),
2463 theASurf2->FirstVParameter()};
2464 const Standard_Real aUppBorder[4] = {theASurf1->LastUParameter(),
2465 theASurf1->LastVParameter(),
2466 theASurf2->LastUParameter(),
2467 theASurf2->LastVParameter()};
2469 // Tune solution tolerance according with object size.
2470 const Standard_Real aRes1 = Max(Precision::PConfusion() / theASurf1->UResolution(1.0),
2471 Precision::PConfusion() / theASurf1->VResolution(1.0));
2472 const Standard_Real aRes2 = Max(Precision::PConfusion() / theASurf2->UResolution(1.0),
2473 Precision::PConfusion() / theASurf2->VResolution(1.0));
2474 const Standard_Real a3DTol = Max(aRes1, aRes2);
2475 const Standard_Real aTol = Max(Precision::Confusion(), a3DTol);
2477 // u1, v1, u2, v2 order is used.
2478 TColStd_Array1OfReal aPnt(1,4);
2479 aPnt(1) = theU1; aPnt(2) = theV1; aPnt(3) = theU2; aPnt(4) = theV2;
2480 TColStd_Array1OfReal aSingularPnt(aPnt);
2482 Standard_Integer aNbIter = 20;
2483 Standard_Boolean aStatus = Standard_False;
2487 aStatus = DistanceMinimizeByGradient(theASurf1, theASurf2, aPnt);
2488 if (aStatus && !AdjustToDomain(4, &aPnt(1), &aLowBorder[0], &aUppBorder[0]))
2491 aStatus = DistanceMinimizeByExtrema(theASurf1, theASurf2->Value(aPnt(3), aPnt(4)),
2493 if (aStatus && !AdjustToDomain(2, &aPnt(1), &aLowBorder[0], &aUppBorder[0]))
2496 aStatus = DistanceMinimizeByExtrema(theASurf2, theASurf1->Value(aPnt(1), aPnt(2)),
2498 if (aStatus && !AdjustToDomain(2, &aPnt(3), &aLowBorder[2], &aUppBorder[2]))
2501 while(!aStatus && (aNbIter > 0));
2503 // Handle singular points.
2504 Standard_Boolean aSingularStatus = HandleSingleSingularPoint(theASurf1,
2508 if (aSingularStatus)
2509 aPnt = aSingularPnt;
2511 if (!aStatus && !aSingularStatus)
2516 gp_Pnt aP1 = theASurf1->Value(aPnt(1), aPnt(2));
2517 gp_Pnt aP2 = theASurf2->Value(aPnt(3), aPnt(4));
2518 const gp_Pnt aPInt(0.5*(aP1.XYZ() + aP2.XYZ()));
2520 const Standard_Real aSQDist = aPInt.SquareDistance(aP1);
2521 if (aSQDist > aTol * aTol)
2526 //Found point is true intersection point
2527 IntSurf_PntOn2S anIP;
2528 anIP.SetValue(aPInt, aPnt(1), aPnt(2), aPnt(3), aPnt(4));
2530 //The main idea of checks below is to define if insertion of
2531 //addition point (on the boundary) does not lead to invalid
2532 //intersection curve (e.g. having a loop).
2534 //Loops are detected with rotation angle of the Walking-line (WL).
2535 //If there is hairpin bend then insertion is forbidden.
2537 //There are at least two possible problems:
2538 // 1. There are some cases when two neighbor points of the WL
2539 // are almost coincident (the distance between them is less
2540 // than Precision::Confusion). It is impossible to define
2541 // rotation angle in these cases. Therefore, points with
2542 // "good" distances should be selected.
2544 // 2. Intersection point on the surface boundary has highest
2545 // priority in compare with other "middle" points. Therefore,
2546 // if insertion of new point will result in a bend then some
2547 // "middle" points should be deleted in order to provide
2548 // correct insertion.
2550 //Problem test cases:
2551 // test bugs modalg_5 bug24585_1
2552 // test boolean bcut_complex G7
2553 // test bugs moddata_2 bug469
2557 while (line->NbPoints() > 1)
2559 const Standard_Integer aNbPnts = line->NbPoints();
2561 Standard_Integer aPInd = 1;
2562 for (; aPInd <= aNbPnts; aPInd++)
2564 aP1.SetXYZ(line->Value(aPInd).Value().XYZ());
2565 if (aP1.SquareDistance(aPInt) > Precision::SquareConfusion())
2569 //else if (aPInd == 1)
2571 // // After insertion, we will obtain
2572 // // two coincident points in the line.
2573 // // Therefore, insertion is forbidden.
2578 for (++aPInd; aPInd <= aNbPnts; aPInd++)
2580 aP2.SetXYZ(line->Value(aPInd).Value().XYZ());
2581 if (aP1.SquareDistance(aP2) > Precision::SquareConfusion())
2585 if (aPInd > aNbPnts)
2590 const gp_XYZ aDir01(aP1.XYZ() - aPInt.XYZ());
2591 const gp_XYZ aDir12(aP2.XYZ() - aP1.XYZ());
2593 if (aDir01.Dot(aDir12) > 0.0)
2601 aP1.SetXYZ(line->Value(1).Value().XYZ());
2602 if (aP1.SquareDistance(aPInt) <= Precision::SquareConfusion())
2607 line->InsertBefore(1, anIP);
2608 isOK = Standard_True;
2612 while (line->NbPoints() > 1)
2614 const Standard_Integer aNbPnts = line->NbPoints();
2616 gp_Pnt aPPrev, aPCurr;
2617 Standard_Integer aPInd = aNbPnts;
2618 for (; aPInd > 0; aPInd--)
2620 aPCurr.SetXYZ(line->Value(aPInd).Value().XYZ());
2621 if (aPCurr.SquareDistance(aPInt) > Precision::SquareConfusion())
2625 //else if (aPInd == aNbPnts)
2627 // // After insertion, we will obtain
2628 // // two coincident points in the line.
2629 // // Therefore, insertion is forbidden.
2634 for (--aPInd; aPInd > 0; aPInd--)
2636 aPPrev.SetXYZ(line->Value(aPInd).Value().XYZ());
2637 if (aPCurr.SquareDistance(aPPrev) > Precision::SquareConfusion())
2646 const gp_XYZ aDirPC(aPCurr.XYZ() - aPPrev.XYZ());
2647 const gp_XYZ aDirCN(aPInt.XYZ() - aPCurr.XYZ());
2649 if (aDirPC.Dot(aDirCN) > 0.0)
2654 RemoveAPoint(aNbPnts);
2657 Standard_Integer aNbPnts = line->NbPoints();
2658 aP1.SetXYZ(line->Value(aNbPnts).Value().XYZ());
2659 if (aP1.SquareDistance(aPInt) <= Precision::SquareConfusion())
2661 RemoveAPoint(aNbPnts);
2665 isOK = Standard_True;
2671 //=======================================================================
2672 //function : PutToBoundary
2674 //=======================================================================
2675 Standard_Boolean IntWalk_PWalking::
2676 PutToBoundary(const Handle(Adaptor3d_Surface)& theASurf1,
2677 const Handle(Adaptor3d_Surface)& theASurf2)
2679 const Standard_Real aTolMin = Precision::Confusion();
2681 Standard_Boolean hasBeenAdded = Standard_False;
2683 const Standard_Real aU1bFirst = theASurf1->FirstUParameter();
2684 const Standard_Real aU1bLast = theASurf1->LastUParameter();
2685 const Standard_Real aU2bFirst = theASurf2->FirstUParameter();
2686 const Standard_Real aU2bLast = theASurf2->LastUParameter();
2687 const Standard_Real aV1bFirst = theASurf1->FirstVParameter();
2688 const Standard_Real aV1bLast = theASurf1->LastVParameter();
2689 const Standard_Real aV2bFirst = theASurf2->FirstVParameter();
2690 const Standard_Real aV2bLast = theASurf2->LastVParameter();
2692 Standard_Real aTol = 1.0;
2693 aTol = Min(aTol, aU1bLast - aU1bFirst);
2694 aTol = Min(aTol, aU2bLast - aU2bFirst);
2695 aTol = Min(aTol, aV1bLast - aV1bFirst);
2696 aTol = Min(aTol, aV2bLast - aV2bFirst)*1.0e-3;
2698 if(aTol <= 2.0*aTolMin)
2699 return hasBeenAdded;
2701 Standard_Boolean isNeedAdding = Standard_False;
2702 Standard_Boolean isU1parallel = Standard_False, isV1parallel = Standard_False;
2703 Standard_Boolean isU2parallel = Standard_False, isV2parallel = Standard_False;
2704 IsParallel(line, Standard_True, aTol, isU1parallel, isV1parallel);
2705 IsParallel(line, Standard_False, aTol, isU2parallel, isV2parallel);
2707 Standard_Real u1, v1, u2, v2;
2708 line->Value(1).Parameters(u1, v1, u2, v2);
2709 Standard_Real aDelta = 0.0;
2713 aDelta = u1 - aU1bFirst;
2714 if((aTolMin < aDelta) && (aDelta < aTol))
2717 isNeedAdding = Standard_True;
2721 aDelta = aU1bLast - u1;
2722 if((aTolMin < aDelta) && (aDelta < aTol))
2725 isNeedAdding = Standard_True;
2732 aDelta = u2 - aU2bFirst;
2733 if((aTolMin < aDelta) && (aDelta < aTol))
2736 isNeedAdding = Standard_True;
2740 aDelta = aU2bLast - u2;
2741 if((aTolMin < aDelta) && (aDelta < aTol))
2744 isNeedAdding = Standard_True;
2751 aDelta = v1 - aV1bFirst;
2752 if((aTolMin < aDelta) && (aDelta < aTol))
2755 isNeedAdding = Standard_True;
2759 aDelta = aV1bLast - v1;
2760 if((aTolMin < aDelta) && (aDelta < aTol))
2763 isNeedAdding = Standard_True;
2770 aDelta = v2 - aV2bFirst;
2771 if((aTolMin < aDelta) && (aDelta < aTol))
2774 isNeedAdding = Standard_True;
2778 aDelta = aV2bLast - v2;
2779 if((aTolMin < aDelta) && (aDelta < aTol))
2782 isNeedAdding = Standard_True;
2790 SeekPointOnBoundary(theASurf1, theASurf2, u1,
2791 v1, u2, v2, Standard_True);
2794 const Standard_Integer aNbPnts = line->NbPoints();
2795 isNeedAdding = Standard_False;
2796 line->Value(aNbPnts).Parameters(u1, v1, u2, v2);
2800 aDelta = u1 - aU1bFirst;
2801 if((aTolMin < aDelta) && (aDelta < aTol))
2804 isNeedAdding = Standard_True;
2808 aDelta = aU1bLast - u1;
2809 if((aTolMin < aDelta) && (aDelta < aTol))
2812 isNeedAdding = Standard_True;
2819 aDelta = u2 - aU2bFirst;
2820 if((aTolMin < aDelta) && (aDelta < aTol))
2823 isNeedAdding = Standard_True;
2827 aDelta = aU2bLast - u2;
2828 if((aTolMin < aDelta) && (aDelta < aTol))
2831 isNeedAdding = Standard_True;
2838 aDelta = v1 - aV1bFirst;
2839 if((aTolMin < aDelta) && (aDelta < aTol))
2842 isNeedAdding = Standard_True;
2846 aDelta = aV1bLast - v1;
2847 if((aTolMin < aDelta) && (aDelta < aTol))
2850 isNeedAdding = Standard_True;
2857 aDelta = v2 - aV2bFirst;
2858 if((aTolMin < aDelta) && (aDelta < aTol))
2861 isNeedAdding = Standard_True;
2865 aDelta = aV2bLast - v2;
2866 if((aTolMin < aDelta) && (aDelta < aTol))
2869 isNeedAdding = Standard_True;
2877 SeekPointOnBoundary(theASurf1, theASurf2, u1,
2878 v1, u2, v2, Standard_False);
2881 return hasBeenAdded;
2884 //=======================================================================
2885 //function : SeekAdditionalPoints
2887 //=======================================================================
2888 Standard_Boolean IntWalk_PWalking::
2889 SeekAdditionalPoints( const Handle(Adaptor3d_Surface)& theASurf1,
2890 const Handle(Adaptor3d_Surface)& theASurf2,
2891 const Standard_Integer theMinNbPoints)
2893 const Standard_Real aTol = 1.0e-14;
2894 Standard_Integer aNbPoints = line->NbPoints();
2895 if(aNbPoints > theMinNbPoints)
2896 return Standard_True;
2898 const Standard_Real aU1bFirst = theASurf1->FirstUParameter();
2899 const Standard_Real aU1bLast = theASurf1->LastUParameter();
2900 const Standard_Real aU2bFirst = theASurf2->FirstUParameter();
2901 const Standard_Real aU2bLast = theASurf2->LastUParameter();
2902 const Standard_Real aV1bFirst = theASurf1->FirstVParameter();
2903 const Standard_Real aV1bLast = theASurf1->LastVParameter();
2904 const Standard_Real aV2bFirst = theASurf2->FirstVParameter();
2905 const Standard_Real aV2bLast = theASurf2->LastVParameter();
2908 Standard_Boolean isPrecise = Standard_False;
2910 TColStd_Array1OfReal aPnt(1, 4);
2913 Standard_Integer aNbPointsPrev = 0;
2914 while(aNbPoints < theMinNbPoints && (aNbPoints != aNbPointsPrev))
2916 aNbPointsPrev = aNbPoints;
2917 for(Standard_Integer fp = 1, lp = 2; fp < aNbPoints; fp = lp + 1)
2919 Standard_Real U1f, V1f, U2f, V2f; //first point in 1st and 2nd surafaces
2920 Standard_Real U1l, V1l, U2l, V2l; //last point in 1st and 2nd surafaces
2923 line->Value(fp).Parameters(U1f, V1f, U2f, V2f);
2924 line->Value(lp).Parameters(U1l, V1l, U2l, V2l);
2926 aPnt(1) = 0.5*(U1f + U1l);
2927 if(aPnt(1) < aU1bFirst)
2928 aPnt(1) = aU1bFirst;
2929 if(aPnt(1) > aU1bLast)
2932 aPnt(2) = 0.5*(V1f+V1l);
2933 if(aPnt(2) < aV1bFirst)
2934 aPnt(2) = aV1bFirst;
2935 if(aPnt(2) > aV1bLast)
2938 aPnt(3) = 0.5*(U2f+U2l);
2939 if(aPnt(3) < aU2bFirst)
2940 aPnt(3) = aU2bFirst;
2941 if(aPnt(3) > aU2bLast)
2944 aPnt(4) = 0.5*(V2f+V2l);
2945 if(aPnt(4) < aV2bFirst)
2946 aPnt(4) = aV2bFirst;
2947 if(aPnt(4) > aV2bLast)
2950 Standard_Boolean aStatus = Standard_False;
2951 Standard_Integer aNbIter = 5;
2954 aStatus = DistanceMinimizeByGradient(theASurf1, theASurf2, aPnt);
2960 aStatus = DistanceMinimizeByExtrema(theASurf1, theASurf2->Value(aPnt(3), aPnt(4)), aPnt(1), aPnt(2));
2966 aStatus = DistanceMinimizeByExtrema(theASurf2, theASurf1->Value(aPnt(1), aPnt(2)), aPnt(3), aPnt(4));
2972 while(!aStatus && (--aNbIter > 0));
2976 gp_Pnt aP1 = theASurf1->Value(aPnt(1), aPnt(2)),
2977 aP2 = theASurf2->Value(aPnt(3), aPnt(4));
2978 gp_Pnt aPInt(0.5*(aP1.XYZ() + aP2.XYZ()));
2980 const Standard_Real aSQDist1 = aPInt.SquareDistance(aP1),
2981 aSQDist2 = aPInt.SquareDistance(aP2);
2983 if((aSQDist1 < aTol) && (aSQDist2 < aTol))
2985 IntSurf_PntOn2S anIP;
2986 anIP.SetValue(aPInt, aPnt(1), aPnt(2), aPnt(3), aPnt(4));
2987 line->InsertBefore(lp, anIP);
2989 isPrecise = Standard_True;
2991 if(++aNbPoints >= theMinNbPoints)
3005 void IntWalk_PWalking::
3006 RepartirOuDiviser(Standard_Boolean& DejaReparti,
3007 IntImp_ConstIsoparametric& ChoixIso,
3008 Standard_Boolean& Arrive)
3010 // at the neighborhood of a point, there is a fail of marching
3011 // it is required to divide the steps to try to continue
3012 // if the step is too small if we are on border
3013 // restart in another direction if it was not done, otherwise stop
3016 // Standard_Integer i;
3017 if (Arrive) { //restart in the other direction
3018 if (!DejaReparti ) {
3019 Arrive = Standard_False;
3020 DejaReparti = Standard_True;
3021 previousPoint = line->Value(1);
3022 previoustg = Standard_False;
3023 previousd1 = firstd1;
3024 previousd2 = firstd2;
3026 myTangentIdx = line->NbPoints();
3030 //-- printf("\nIntWalk_PWalking_2.gxx Reverse %3d\n",indextg);
3031 sensCheminement = -1;
3033 tglast = Standard_False;
3034 ChoixIso = choixIsoSav;
3041 Standard_Real u1,v1,u2,v2;
3042 Standard_Real U1,V1,U2,V2;
3043 Standard_Integer nn=line->NbPoints();
3045 line->Value(nn).Parameters(u1,v1,u2,v2);
3046 line->Value(nn-1).Parameters(U1,V1,U2,V2);
3047 pasuv[0]=Abs(u1-U1);
3048 pasuv[1]=Abs(v1-V1);
3049 pasuv[2]=Abs(u2-U2);
3050 pasuv[3]=Abs(v2-V2);
3057 if ( pasuv[0]*0.5 < ResoU1
3058 && pasuv[1]*0.5 < ResoV1
3059 && pasuv[2]*0.5 < ResoU2
3060 && pasuv[3]*0.5 < ResoV2
3063 tglast = Standard_True; // IS IT ENOUGH ????
3067 { //restart in the other direction
3068 DejaReparti = Standard_True;
3069 previousPoint = line->Value(1);
3070 previoustg = Standard_False;
3071 previousd1 = firstd1;
3072 previousd2 = firstd2;
3074 myTangentIdx = line->NbPoints();
3078 //-- printf("\nIntWalk_PWalking_2.gxx Reverse %3d\n",indextg);
3080 sensCheminement = -1;
3082 tglast = Standard_False;
3083 ChoixIso = choixIsoSav;
3091 Standard_Real u1,v1,u2,v2;
3092 Standard_Real U1,V1,U2,V2;
3093 Standard_Integer nn=line->NbPoints();
3095 line->Value(nn).Parameters(u1,v1,u2,v2);
3096 line->Value(nn-1).Parameters(U1,V1,U2,V2);
3097 pasuv[0]=Abs(u1-U1);
3098 pasuv[1]=Abs(v1-V1);
3099 pasuv[2]=Abs(u2-U2);
3100 pasuv[3]=Abs(v2-V2);
3104 else Arrive = Standard_True;
3116 //OCC431(apo): modified ->
3117 static const Standard_Real CosRef2D = Cos(M_PI/9.0), AngRef2D = M_PI/2.0;
3119 static const Standard_Real d = 7.0;
3122 IntWalk_StatusDeflection IntWalk_PWalking::TestDeflection(const IntImp_ConstIsoparametric choixIso,
3123 const IntWalk_StatusDeflection theStatus)
3125 // test if vector is observed by calculating an increase of vector
3126 // or the previous point and its tangent, the new calculated point and its
3127 // tangent; it is possible to find a cube passing by the 2 points and having as a
3128 // derivative the tangents of the intersection
3129 // calculate the point with parameter 0.5 on cube=p1
3130 // calculate the medium point of 2 points of intersection=p2
3131 // if arrow/2<=||p1p2||<= arrow consider that the vector is observed
3132 // otherwise adjust the step depending on the ratio ||p1p2||/vector
3133 // and the previous step
3134 // test if in 2 tangent planes of surfaces there is no too great angle2d
3135 // grand : if yes divide the step
3136 // test if there is no change of side
3139 if(line->NbPoints() ==1 ) {
3140 STATIC_BLOCAGE_SUR_PAS_TROP_GRAND=STATIC_PRECEDENT_INFLEXION=0;
3143 IntWalk_StatusDeflection aStatus = IntWalk_OK;
3144 Standard_Real FlecheCourante , Ratio = 1.0;
3147 const Handle(Adaptor3d_Surface)& Caro1 = myIntersectionOn2S.Function().AuxillarSurface1();
3148 const Handle(Adaptor3d_Surface)& Caro2 = myIntersectionOn2S.Function().AuxillarSurface2();
3150 const IntSurf_PntOn2S& CurrentPoint = myIntersectionOn2S.Point();
3151 //==================================================================================
3152 //========= S t o p o n p o i n t ============
3153 //==================================================================================
3154 if (myIntersectionOn2S.IsTangent()) {
3155 return IntWalk_ArretSurPoint;
3158 const gp_Dir& TgCourante = myIntersectionOn2S.Direction();
3160 const Standard_Real aCosBetweenTangent = TgCourante.Dot(previousd);
3162 //==================================================================================
3163 //========= R i s k o f i n f l e x i o n p o i n t ============
3164 //==================================================================================
3165 if (aCosBetweenTangent < 0) {
3166 //------------------------------------------------------------
3167 //-- Risk of inflexion point : Divide the step by 2
3168 //-- Initialize STATIC_PRECEDENT_INFLEXION so that
3169 //-- at the next call to return Pas_OK if there is no
3170 //-- more risk of the point of inflexion
3171 //------------------------------------------------------------
3177 STATIC_PRECEDENT_INFLEXION+=3;
3178 if (pasuv[0] < ResoU1 && pasuv[1] <ResoV1 && pasuv[2] <ResoU2 && pasuv[3] < ResoV2)
3179 return IntWalk_ArretSurPointPrecedent;
3181 return IntWalk_PasTropGrand;
3184 if(STATIC_PRECEDENT_INFLEXION > 0) {
3185 STATIC_PRECEDENT_INFLEXION -- ;
3190 //==================================================================================
3191 //========= D e t e c t c o n f u s e d P o in t s ===========
3192 //==================================================================================
3194 const Standard_Real aSqDist = previousPoint.Value().
3195 SquareDistance(CurrentPoint.Value());
3198 if (aSqDist < Precision::SquareConfusion()) {
3199 pasInit[0] = Max(pasInit[0], 5.0*ResoU1);
3200 pasInit[1] = Max(pasInit[1], 5.0*ResoV1);
3201 pasInit[2] = Max(pasInit[2], 5.0*ResoU2);
3202 pasInit[3] = Max(pasInit[3], 5.0*ResoV2);
3204 for(Standard_Integer i = 0; i < 4; i++)
3206 pasuv[i] = Max(pasuv[i], Min(1.5*pasuv[i], pasInit[i]));
3208 //Compute local resolution: for OCC26717
3209 if (Abs(pasuv[choixIso] - pasInit[choixIso]) <= Precision::Confusion())
3211 Standard_Real CurU, CurV;
3212 if (choixIso == IntImp_UIsoparametricOnCaro1 ||
3213 choixIso == IntImp_VIsoparametricOnCaro1)
3214 previousPoint.ParametersOnS1(CurU, CurV);
3216 previousPoint.ParametersOnS2(CurU, CurV);
3217 gp_Pnt CurPnt = (choixIso == IntImp_UIsoparametricOnCaro1 ||
3218 choixIso == IntImp_VIsoparametricOnCaro1)?
3219 Adaptor3d_HSurfaceTool::Value(Caro1, CurU, CurV) :
3220 Adaptor3d_HSurfaceTool::Value(Caro2, CurU, CurV);
3224 case IntImp_UIsoparametricOnCaro1:
3226 Adaptor3d_HSurfaceTool::Value(Caro1,
3227 CurU + sensCheminement*pasuv[0],
3230 case IntImp_VIsoparametricOnCaro1:
3232 Adaptor3d_HSurfaceTool::Value(Caro1,
3234 CurV + sensCheminement*pasuv[1]);
3236 case IntImp_UIsoparametricOnCaro2:
3238 Adaptor3d_HSurfaceTool::Value(Caro2,
3239 CurU + sensCheminement*pasuv[2],
3242 case IntImp_VIsoparametricOnCaro2:
3244 Adaptor3d_HSurfaceTool::Value(Caro2,
3246 CurV + sensCheminement*pasuv[3]);
3250 Standard_Real RefDist = CurPnt.Distance(OffsetPnt);
3251 Standard_Real LocalResol = 0.;
3252 if (RefDist > gp::Resolution())
3253 LocalResol = pasuv[choixIso] * tolconf / RefDist;
3254 if (pasuv[choixIso] < 2*LocalResol)
3255 pasuv[choixIso] = pasInit[choixIso] = 2*LocalResol;
3257 ////////////////////////////////////////
3258 aStatus = IntWalk_PointConfondu;
3261 //==================================================================================
3262 Standard_Real Up1,Vp1,Uc1,Vc1,Du1,Dv1,AbsDu1,AbsDu2,AbsDv1,AbsDv2;
3263 Standard_Real Up2,Vp2,Uc2,Vc2,Du2,Dv2;
3265 previousPoint.Parameters(Up1,Vp1,Up2,Vp2);
3266 CurrentPoint.Parameters(Uc1,Vc1,Uc2,Vc2);
3268 Du1 = Uc1 - Up1; Dv1 = Vc1 - Vp1;
3269 Du2 = Uc2 - Up2; Dv2 = Vc2 - Vp2;
3275 //=================================================================================
3276 //==== S t e p o f p r o g r e s s i o n (between previous and Current) =======
3277 //=================================================================================
3278 if ( AbsDu1 < ResoU1 && AbsDv1 < ResoV1
3279 && AbsDu2 < ResoU2 && AbsDv2 < ResoV2) {
3280 pasuv[0] = ResoU1; pasuv[1] = ResoV1; pasuv[2] = ResoU2; pasuv[3] = ResoV2;
3281 return(IntWalk_ArretSurPointPrecedent);
3283 //==================================================================================
3285 Standard_Real tolArea = 100.0;
3286 if (ResoU1 < Precision::PConfusion() ||
3287 ResoV1 < Precision::PConfusion() ||
3288 ResoU2 < Precision::PConfusion() ||
3289 ResoV2 < Precision::PConfusion() )
3290 tolArea = tolArea*2.0;
3292 Standard_Real Cosi1, CosRef1, Ang1, AngRef1, ResoUV1, Duv1, d1, tolCoeff1;
3293 Standard_Real Cosi2, CosRef2, Ang2, AngRef2, ResoUV2, Duv2, d2, tolCoeff2;
3294 Cosi1 = Du1*previousd1.X() + Dv1*previousd1.Y();
3295 Cosi2 = Du2*previousd2.X() + Dv2*previousd2.Y();
3296 Duv1 = Du1*Du1 + Dv1*Dv1;
3297 Duv2 = Du2*Du2 + Dv2*Dv2;
3298 ResoUV1 = ResoU1*ResoU1 + ResoV1*ResoV1;
3299 ResoUV2 = ResoU2*ResoU2 + ResoV2*ResoV2;
3301 //modified by NIZNHY-PKV Wed Nov 13 12:25:44 2002 f
3303 Standard_Real aMinDiv2=Precision::Confusion();
3304 aMinDiv2=aMinDiv2*aMinDiv2;
3307 if (Duv1>aMinDiv2) {
3308 d1 = Abs(ResoUV1/Duv1);
3309 d1 = Min(Sqrt(d1)*tolArea, d);
3311 //d1 = Abs(ResoUV1/Duv1);
3312 //d1 = Min(Sqrt(d1)*tolArea,d);
3313 //modified by NIZNHY-PKV Wed Nov 13 12:34:30 2002 t
3314 tolCoeff1 = Exp(d1);
3316 //modified by NIZNHY-PKV Wed Nov 13 12:34:43 2002 f
3318 if (Duv2>aMinDiv2) {
3319 d2 = Abs(ResoUV2/Duv2);
3320 d2 = Min(Sqrt(d2)*tolArea,d);
3322 //d2 = Abs(ResoUV2/Duv2);
3323 //d2 = Min(Sqrt(d2)*tolArea,d);
3324 //modified by NIZNHY-PKV Wed Nov 13 12:34:53 2002 t
3325 tolCoeff2 = Exp(d2);
3326 CosRef1 = CosRef2D/tolCoeff1;
3327 CosRef2 = CosRef2D/tolCoeff2;
3329 //==================================================================================
3330 //== The points are not confused : ==
3331 //== D e t e c t t h e S t o p a t p r e v i o u s p o i n t ==
3332 //== N o t T o o G r e a t (angle in space UV) ==
3333 //== C h a n g e o f s i d e ==
3334 //==================================================================================
3335 if (aStatus != IntWalk_PointConfondu) {
3336 if(Cosi1*Cosi1 < CosRef1*Duv1 || Cosi2*Cosi2 < CosRef2*Duv2) {
3337 pasuv[0]*=0.5; pasuv[1]*=0.5; pasuv[2]*=0.5; pasuv[3]*=0.5;
3338 if (pasuv[0]<ResoU1 && pasuv[1]<ResoV1 && pasuv[2]<ResoU2 && pasuv[3]<ResoV2) {
3339 return(IntWalk_ArretSurPointPrecedent);
3342 pasuv[0]*=0.5; pasuv[1]*=0.5; pasuv[2]*=0.5; pasuv[3]*=0.5;
3343 return(IntWalk_PasTropGrand);
3346 const gp_Dir2d& Tg2dcourante1 = myIntersectionOn2S.DirectionOnS1();
3347 const gp_Dir2d& Tg2dcourante2 = myIntersectionOn2S.DirectionOnS2();
3348 Cosi1 = Du1*Tg2dcourante1.X() + Dv1*Tg2dcourante1.Y();
3349 Cosi2 = Du2*Tg2dcourante2.X() + Dv2*Tg2dcourante2.Y();
3350 Ang1 = Abs(previousd1.Angle(Tg2dcourante1));
3351 Ang2 = Abs(previousd2.Angle(Tg2dcourante2));
3352 AngRef1 = AngRef2D*tolCoeff1;
3353 AngRef2 = AngRef2D*tolCoeff2;
3354 //-------------------------------------------------------
3355 //-- Test : Angle too great in space UV -----
3356 //-- Change of side -----
3357 //-------------------------------------------------------
3358 if(Cosi1*Cosi1 < CosRef1*Duv1 || Cosi2*Cosi2 < CosRef2*Duv2 || Ang1 > AngRef1 || Ang2 > AngRef2) {
3359 pasuv[0]*=0.5; pasuv[1]*=0.5; pasuv[2]*=0.5; pasuv[3]*=0.5;
3360 if (pasuv[0]<ResoU1 && pasuv[1]<ResoV1 && pasuv[2]<ResoU2 && pasuv[3]<ResoV2)
3361 return(IntWalk_ArretSurPoint);
3363 return(IntWalk_PasTropGrand);
3367 //==================================================================================
3368 //== D e t e c t i o n o f : Step Too Small
3370 //==================================================================================
3372 //---------------------------------------
3373 //-- Estimate of the vector --
3374 //---------------------------------------
3376 Sqrt(Abs((previousd.XYZ()-TgCourante.XYZ()).SquareModulus()*aSqDist))/8.;
3378 if ( FlecheCourante<= fleche*0.5) { //-- Current step too small
3379 if(FlecheCourante>1e-16) {
3380 Ratio = 0.5*(fleche/FlecheCourante);
3385 Standard_Real pasSu1 = pasuv[0];
3386 Standard_Real pasSv1 = pasuv[1];
3387 Standard_Real pasSu2 = pasuv[2];
3388 Standard_Real pasSv2 = pasuv[3];
3391 //-- a point at U+DeltaU is required, ....
3392 //-- return a point at U + Epsilon
3393 //-- Epsilon << DeltaU.
3395 if(pasuv[0]< AbsDu1) pasuv[0] = AbsDu1;
3396 if(pasuv[1]< AbsDv1) pasuv[1] = AbsDv1;
3397 if(pasuv[2]< AbsDu2) pasuv[2] = AbsDu2;
3398 if(pasuv[3]< AbsDv2) pasuv[3] = AbsDv2;
3400 if(pasuv[0]<ResoU1) pasuv[0]=ResoU1;
3401 if(pasuv[1]<ResoV1) pasuv[1]=ResoV1;
3402 if(pasuv[2]<ResoU2) pasuv[2]=ResoU2;
3403 if(pasuv[3]<ResoV2) pasuv[3]=ResoV2;
3404 //-- if(Ratio>10.0 ) { Ratio=10.0; }
3405 Standard_Real R1,R = pasInit[0]/pasuv[0];
3406 R1= pasInit[1]/pasuv[1]; if(R1<R) R=R1;
3407 R1= pasInit[2]/pasuv[2]; if(R1<R) R=R1;
3408 R1= pasInit[3]/pasuv[3]; if(R1<R) R=R1;
3409 if(Ratio > R) Ratio=R;
3410 pasuv[0] = Min(Ratio*pasuv[0],pasInit[0]);
3411 pasuv[1] = Min(Ratio*pasuv[1],pasInit[1]);
3412 pasuv[2] = Min(Ratio*pasuv[2],pasInit[2]);
3413 pasuv[3] = Min(Ratio*pasuv[3],pasInit[3]);
3414 if (pasuv[0] != pasSu1 || pasuv[2] != pasSu2||
3415 pasuv[1] != pasSv1 || pasuv[3] != pasSv2) {
3416 if(++STATIC_BLOCAGE_SUR_PAS_TROP_GRAND > 5) {
3417 STATIC_BLOCAGE_SUR_PAS_TROP_GRAND = 0;
3418 return IntWalk_PasTropGrand;
3421 if(aStatus == IntWalk_OK) {
3422 STATIC_BLOCAGE_SUR_PAS_TROP_GRAND=0;
3423 //-- Try to increase the step
3427 else { //-- CurrentVector > vector*0.5
3428 if (FlecheCourante > fleche) { //-- Current step too Great
3429 Ratio = fleche/FlecheCourante;
3430 pasuv[0] = Ratio*pasuv[0];
3431 pasuv[1] = Ratio*pasuv[1];
3432 pasuv[2] = Ratio*pasuv[2];
3433 pasuv[3] = Ratio*pasuv[3];
3434 //if(++STATIC_BLOCAGE_SUR_PAS_TROP_GRAND > 5) {
3435 // STATIC_BLOCAGE_SUR_PAS_TROP_GRAND = 0;
3436 return IntWalk_PasTropGrand;
3439 else { //-- vector/2 < CurrentVector <= vector
3440 Ratio = 0.75 * (fleche / FlecheCourante);
3444 if(aStatus != IntWalk_PointConfondu)
3446 //Here, aCosBetweenTangent >= 0.0 definitely.
3449 Brief algorithm description.
3450 We have two (not-coincindent) intersection point (P1 and P2). In every point,
3451 vector of tangent (to the intersection curve) is known (vectors T1 and T2).
3452 Basing on these data, we create osculating circle.
3454 * - arc of osculating circle
3461 Let me draw your attention to the following facts:
3462 1. Vectors T1 and T2 direct FROM (not TO) points P1 and P2. Therefore,
3463 one of previously computed vector should be reversed.
3465 In this case, the absolute (!) value of the deflection between the arc of
3466 the osculating circle and the P1P2 segment can be computed as follows:
3467 e = d*(1-sin(B/2))/(2*cos(B/2)), (1)
3468 where d is the length of P1P2 segment, B is the angle between vectors T1 and T2.
3476 Later, the normal state of algorithm work is (as we apply)
3477 tolconf/2 <= e <= tolconf.
3478 In this case, we keep previous step.
3480 If e < tolconf/2 then the local curvature of the intersection curve is small.
3481 As result, the step should be increased.
3483 If e > tolconf then the step is too big. Therefore, we should decrease one.
3485 Condition (1) is equivalent to
3486 sin(B/2) = 1 - 2/(1+(d/(2*e))^2) = Fs(e),
3487 cos(B) = 1 - 2*Fs(e)^2 = Fd(e),
3488 where Fs(e)and Fd(e) are some function with parameter "deflection".
3490 Let mean that Fs(e) is decreasing function. Fd(e) is increasing function,
3491 in the range, where Fs(e) > 0.0 (i.e. when e < d/2).
3493 Now, let substitute required deflection (tolconf or tolconf/2) to z. Then
3494 it is necessary to check if e < z or if e > z.
3496 In this case, it is enough to comapare Fs(e) and Fs(z).
3497 At that Fs(e) > 0 because sin(B/2) > 0 always.
3499 Therefore, if Fs(z) < 0.0 then Fs(e) > Fs(z) ==> e < z definitely.
3500 If Fs(z) > 0.0 then we can compare Fs(z)^2 and Fs(e)^2 or, in substance,
3501 values Fd(e) and Fd(z). If Fd(e) > Fd(z) then e > z and vice versa.
3504 //Fd(e) is already known (Fd(e) == -aCosBetweenTangent)
3506 const Standard_Real anInvSqAbsArcDeflMax = 0.25*aSqDist/(tolconf*tolconf);
3507 const Standard_Real aSinB2Max = 1.0 - 2.0/(1.0 + anInvSqAbsArcDeflMax);
3509 if(aSinB2Max >= 0.0 && (aCosBetweenTangent <= 2.0 * aSinB2Max * aSinB2Max - 1.0))
3510 {//Real deflection is greater or equal than tolconf
3511 aStatus = IntWalk_PasTropGrand;
3514 {//Real deflection is less than tolconf
3515 const Standard_Real anInvSqAbsArcDeflMin = 4.0*anInvSqAbsArcDeflMax;
3516 const Standard_Real aSinB2Min = 1.0 - 2.0/(1.0 + anInvSqAbsArcDeflMin);
3518 if (theStatus != IntWalk_PasTropGrand &&
3519 ((aSinB2Min < 0.0) || (aCosBetweenTangent >= 2.0 * aSinB2Min * aSinB2Min - 1.0)))
3520 {//Real deflection is less than tolconf/2.0
3521 aStatus = IntWalk_StepTooSmall;
3525 if(aStatus == IntWalk_PasTropGrand)
3527 pasuv[0]*=0.5; pasuv[1]*=0.5; pasuv[2]*=0.5; pasuv[3]*=0.5;
3531 if(aStatus == IntWalk_StepTooSmall)
3533 pasuv[0] = Max(pasuv[0], AbsDu1);
3534 pasuv[1] = Max(pasuv[1], AbsDv1);
3535 pasuv[2] = Max(pasuv[2], AbsDu2);
3536 pasuv[3] = Max(pasuv[3], AbsDv2);
3538 pasInit[0] = Max(pasInit[0], AbsDu1);
3539 pasInit[1] = Max(pasInit[1], AbsDv1);
3540 pasInit[2] = Max(pasInit[2], AbsDu2);
3541 pasInit[3] = Max(pasInit[3], AbsDv2);
3547 pasuv[0] = Max(myStepMin[0],Min(Min(Ratio*AbsDu1,pasuv[0]),pasInit[0]));
3548 pasuv[1] = Max(myStepMin[1],Min(Min(Ratio*AbsDv1,pasuv[1]),pasInit[1]));
3549 pasuv[2] = Max(myStepMin[2],Min(Min(Ratio*AbsDu2,pasuv[2]),pasInit[2]));
3550 pasuv[3] = Max(myStepMin[3],Min(Min(Ratio*AbsDv2,pasuv[3]),pasInit[3]));
3552 if(aStatus == IntWalk_OK) STATIC_BLOCAGE_SUR_PAS_TROP_GRAND=0;
3556 Standard_Boolean IntWalk_PWalking::
3557 TestArret(const Standard_Boolean DejaReparti,
3558 TColStd_Array1OfReal& Param,
3559 IntImp_ConstIsoparametric& ChoixIso)
3562 // test if the point of intersection set by these parameters remains in the
3563 // natural domain of each square.
3564 // if the point outpasses reframe to find the best iso (border)
3565 // that intersects easiest the other square
3566 // otherwise test if closed line is present
3569 Standard_Real Uvd[4],Uvf[4],Epsuv[4],Duv[4],Uvp[4],dv,dv2,ParC[4];
3570 Standard_Real DPc,DPb;
3571 Standard_Integer i = 0, k = 0;
3576 previousPoint.Parameters(Uvp[0],Uvp[1],Uvp[2],Uvp[3]);
3578 Standard_Real SolParam[4];
3579 myIntersectionOn2S.Point().Parameters(SolParam[0],SolParam[1],SolParam[2],SolParam[3]);
3581 Standard_Boolean Trouve = Standard_False;
3583 Uvd[0]=Um1; Uvf[0]=UM1; Uvd[1]=Vm1; Uvf[1]=VM1;
3584 Uvd[2]=Um2; Uvf[2]=UM2; Uvd[3]=Vm2; Uvf[3]=VM2;
3586 Standard_Integer im1;
3587 for ( i = 1,im1 = 0;i<=4;i++,im1++) {
3594 if (Param(i) < (Uvd[im1]-Epsuv[im1]) ||
3595 SolParam[im1] < (Uvd[im1]-Epsuv[im1])) //-- Current ----- Bound Inf ----- Previous
3597 Trouve = Standard_True; //--
3598 DPc = Uvp[im1]-Param(i); //-- Previous - Current
3599 DPb = Uvp[im1]-Uvd[im1]; //-- Previous - Bound Inf
3600 ParC[im1] = Uvd[im1]; //-- ParamCorrige
3601 dv = Param(k)-Uvp[k-1]; //-- Current - Previous (other Direction)
3603 if(dv2>RealEpsilon()) { //-- Progress at the other Direction ?
3604 Duv[im1] = DPc*DPb + dv2;
3605 Duv[im1] = Duv[im1]*Duv[im1]/(DPc*DPc+dv2)/(DPb*DPb+dv2);
3608 Duv[im1]=-1.0; //-- If no progress, do not change
3609 } //-- the choice of iso
3611 else if (Param(i) > (Uvf[im1] + Epsuv[im1]) ||
3612 SolParam[im1] > (Uvf[im1] + Epsuv[im1]))//-- Previous ----- Bound Sup ----- Current
3614 Trouve = Standard_True; //--
3615 DPc = Param(i)-Uvp[im1]; //-- Current - Previous
3616 DPb = Uvf[im1]-Uvp[im1]; //-- Bound Sup - Previous
3617 ParC[im1] = Uvf[im1]; //-- Param Corrige
3618 dv = Param(k)-Uvp[k-1]; //-- Current - Previous (other Direction)
3620 if(dv2>RealEpsilon()) { //-- Progress in other Direction ?
3621 Duv[im1] = DPc*DPb + dv2;
3622 Duv[im1] = Duv[im1]*Duv[im1]/(DPc*DPc+dv2)/(DPb*DPb+dv2);
3625 Duv[im1]=-1.0; //-- If no progress, do not change
3626 } //-- the choice of iso
3635 //--------------------------------------------------
3636 //-- One of Parameters u1,v1,u2,v2 is outside of --
3637 //-- the natural limits. --
3638 //-- Find the best direction of --
3639 //-- progress and reframe the parameters. --
3640 //--------------------------------------------------
3641 Standard_Real ddv = -1.0;
3643 for (i=0;i<=3;i++) {
3644 Param(i+1) = ParC[i];
3651 ChoixIso = ChoixRef(k);
3654 if((ParC[0]<=Uvd[0]+Epsuv[0]) || (ParC[0]>=Uvf[0]-Epsuv[0])) {
3655 ChoixIso = IntImp_UIsoparametricOnCaro1;
3657 else if((ParC[1]<=Uvd[1]+Epsuv[1]) || (ParC[1]>=Uvf[1]-Epsuv[1])) {
3658 ChoixIso = IntImp_VIsoparametricOnCaro1;
3660 else if((ParC[2]<=Uvd[2]+Epsuv[2]) || (ParC[2]>=Uvf[2]-Epsuv[2])) {
3661 ChoixIso = IntImp_UIsoparametricOnCaro2;
3663 else if((ParC[3]<=Uvd[3]+Epsuv[3]) || (ParC[3]>=Uvf[3]-Epsuv[3])) {
3664 ChoixIso = IntImp_VIsoparametricOnCaro2;
3667 close = Standard_False;
3668 return Standard_True;
3672 if (!DejaReparti) { // find if line closed
3675 const IntSurf_PntOn2S& POn2S1=line->Value(1);
3677 POn2S1.ParametersOnS1(u,v);
3678 gp_Pnt2d P1uvS1(u,v);
3679 previousPoint.ParametersOnS1(u,v);
3680 gp_Pnt2d PrevuvS1(u,v);
3681 myIntersectionOn2S.Point().ParametersOnS1(u,v);
3682 gp_Pnt2d myIntersuvS1(u,v);
3683 Standard_Boolean close2dS1 = (P1uvS1.XY()-PrevuvS1.XY())*
3684 (P1uvS1.XY()-myIntersuvS1.XY()) < 0.0;
3686 POn2S1.ParametersOnS2(u,v);
3687 gp_Pnt2d P1uvS2(u,v);
3688 previousPoint.ParametersOnS2(u,v);
3689 gp_Pnt2d PrevuvS2(u,v);
3690 myIntersectionOn2S.Point().ParametersOnS2(u,v);
3691 gp_Pnt2d myIntersuvS2(u,v);
3692 Standard_Boolean close2dS2 = (P1uvS2.XY()-PrevuvS2.XY())*
3693 (P1uvS2.XY()-myIntersuvS2.XY()) < 0.0;
3695 close = close2dS1 && close2dS2;
3698 else return Standard_False;