0026151: Wrong result obtained by intersection algorithm.
[occt.git] / src / IntWalk / IntWalk_PWalking.cxx
1 // Copyright (c) 1995-1999 Matra Datavision
2 // Copyright (c) 1999-2014 OPEN CASCADE SAS
3 //
4 // This file is part of Open CASCADE Technology software library.
5 //
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.
11 //
12 // Alternatively, this file may be used under the terms of Open CASCADE
13 // commercial license or contractual agreement.
14
15 #include <IntWalk_PWalking.ixx>
16
17 #include <IntWalk_StatusDeflection.hxx>
18
19 #include <TColgp_Array1OfPnt.hxx>
20 #include <TColStd_Array1OfReal.hxx>
21
22 #include <IntImp_ComputeTangence.hxx>
23
24 #include <Adaptor3d_HSurface.hxx>
25 #include <Adaptor3d_HSurfaceTool.hxx>
26
27 #include <Precision.hxx>
28
29 #include <math_FunctionSetRoot.hxx>
30 #include <Geom_Surface.hxx>
31
32 #include <Standard_Failure.hxx>
33 #include <gp_Pnt2d.hxx>
34
35 #include <Extrema_GenLocateExtPS.hxx>
36
37 //==================================================================================
38 // function : IntWalk_PWalking::IntWalk_PWalking
39 // purpose  :
40 // estimate of max step : To avoid abrupt changes 
41 // during change of isos 
42 //==================================================================================
43 void ComputePasInit(Standard_Real *pasuv,
44                     Standard_Real Um1,Standard_Real UM1,
45                     Standard_Real Vm1,Standard_Real VM1,
46                     Standard_Real Um2,Standard_Real UM2,
47                     Standard_Real Vm2,Standard_Real VM2,
48                     Standard_Real _Um1,Standard_Real _UM1,
49                     Standard_Real _Vm1,Standard_Real _VM1,
50                     Standard_Real _Um2,Standard_Real _UM2,
51                     Standard_Real _Vm2,Standard_Real _VM2,
52                     const Handle(Adaptor3d_HSurface)& ,
53                     const Handle(Adaptor3d_HSurface)& ,
54                     const Standard_Real Increment) 
55
56   Standard_Real du1=Abs(UM1-Um1);
57   Standard_Real dv1=Abs(VM1-Vm1);
58   Standard_Real du2=Abs(UM2-Um2);
59   Standard_Real dv2=Abs(VM2-Vm2);
60
61   Standard_Real _du1=Abs(_UM1-_Um1);
62   Standard_Real _dv1=Abs(_VM1-_Vm1);
63   Standard_Real _du2=Abs(_UM2-_Um2);
64   Standard_Real _dv2=Abs(_VM2-_Vm2);
65
66   //-- limit the reduction of uv box estimate to 0.01 natural box
67   //--  du1 : On box of Inter
68   //-- _du1 : On parametric space
69   if(_du1<1e50 && du1<0.01*_du1) du1=0.01*_du1;
70   if(_dv1<1e50 && dv1<0.01*_dv1) dv1=0.01*_dv1;
71   if(_du2<1e50 && du2<0.01*_du2) du2=0.01*_du2;
72   if(_dv2<1e50 && dv2<0.01*_dv2) dv2=0.01*_dv2;
73
74   pasuv[0]=Increment*du1;
75   pasuv[1]=Increment*dv1;
76   pasuv[2]=Increment*du2;
77   pasuv[3]=Increment*dv2;
78 }
79
80 //=======================================================================
81 //function : IsParallel
82 //purpose  : Checks if theLine is parallel of some boundary of given
83 //            surface (it is determined by theCheckSurf1 flag). 
84 //            Parallelism assumes small oscillations (swing is less or 
85 //            equal than theToler).
86 //            Small lines (if first and last parameters in the Surface 
87 //            are almost equal) are classified as parallel (as same as 
88 //            any point can be considered as parallel of any line).
89 //=======================================================================
90 static void IsParallel(const Handle(IntSurf_LineOn2S)& theLine,
91                        const Standard_Boolean theCheckSurf1,
92                        const Standard_Real theToler,
93                        Standard_Boolean& theIsUparallel,
94                        Standard_Boolean& theIsVparallel)
95 {
96   const Standard_Integer aNbPointsMAX = 23;
97
98   theIsUparallel = theIsVparallel = Standard_True;
99
100   Standard_Integer aNbPoints = theLine->NbPoints();
101   if(aNbPoints > aNbPointsMAX)
102   {
103     aNbPoints = aNbPointsMAX;
104   }
105   else if(aNbPoints < 3)
106   {
107     //Here we cannot estimate parallelism.
108     //Do all same as for small lines 
109     return;
110   }
111
112   Standard_Real aStep = IntToReal(theLine->NbPoints()) / aNbPoints;
113   Standard_Real aNPoint = 1.0;
114
115   Standard_Real aUmin = RealLast(), aUmax = RealFirst(), aVmin = RealLast(), aVmax = RealFirst();
116   for(Standard_Integer aNum = 1; aNum <= aNbPoints; aNum++, aNPoint += aStep)
117   {
118     if(aNPoint > aNbPoints)
119     {
120       aNPoint = aNbPoints;
121     }
122
123     Standard_Real u, v;
124     if(theCheckSurf1)
125       theLine->Value(RealToInt(aNPoint)).ParametersOnS1(u, v);
126     else
127       theLine->Value(RealToInt(aNPoint)).ParametersOnS2(u, v);
128
129     if(u < aUmin)
130       aUmin = u;
131
132     if(u > aUmax)
133       aUmax = u;
134
135     if(v < aVmin)
136       aVmin = v;
137
138     if(v > aVmax)
139       aVmax = v;
140   }
141
142   theIsVparallel = ((aUmax - aUmin) < theToler);
143   theIsUparallel = ((aVmax - aVmin) < theToler);
144 }
145
146 //=======================================================================
147 //function : Checking
148 //purpose  : Check, if given point is in surface's boundaries.
149 //            If "yes" then theFactTol = 0.0, else theFactTol is
150 //            equal maximal deviation.
151 //=======================================================================
152 static Standard_Boolean Checking( const Handle(Adaptor3d_HSurface)& theASurf1,
153                                  const Handle(Adaptor3d_HSurface)& theASurf2,
154                                  Standard_Real& theU1,
155                                  Standard_Real& theV1,
156                                  Standard_Real& theU2,
157                                  Standard_Real& theV2,
158                                  Standard_Real& theFactTol)
159 {
160   const Standard_Real aTol = Precision::PConfusion();
161   const Standard_Real aU1bFirst = theASurf1->FirstUParameter();
162   const Standard_Real aU1bLast = theASurf1->LastUParameter();
163   const Standard_Real aU2bFirst = theASurf2->FirstUParameter();
164   const Standard_Real aU2bLast = theASurf2->LastUParameter();
165   const Standard_Real aV1bFirst = theASurf1->FirstVParameter();
166   const Standard_Real aV1bLast = theASurf1->LastVParameter();
167   const Standard_Real aV2bFirst = theASurf2->FirstVParameter();
168   const Standard_Real aV2bLast = theASurf2->LastVParameter();
169
170   Standard_Boolean isOnOrIn = Standard_True;
171   theFactTol = 0.0;
172
173   Standard_Real aDelta = aU1bFirst - theU1;
174   if(aDelta > aTol)
175   {
176     theU1 = aU1bFirst;
177     theFactTol = Max(theFactTol, aDelta);
178     isOnOrIn = Standard_False;
179   }
180
181   aDelta = theU1 - aU1bLast;
182   if(aDelta > aTol)
183   {
184     theU1 = aU1bLast;
185     theFactTol = Max(theFactTol, aDelta);
186     isOnOrIn = Standard_False;
187   }
188
189   aDelta = aV1bFirst - theV1;
190   if(aDelta > aTol)
191   {
192     theV1 = aV1bFirst;
193     theFactTol = Max(theFactTol, aDelta);
194     isOnOrIn = Standard_False;
195   }
196
197   aDelta = theV1 - aV1bLast;
198   if(aDelta > aTol)
199   {
200     theV1 = aV1bLast;
201     theFactTol = Max(theFactTol, aDelta);
202     isOnOrIn = Standard_False;
203   }
204
205   aDelta = aU2bFirst - theU2;
206   if(aDelta > aTol)
207   {
208     theU2 = aU2bFirst;
209     theFactTol = Max(theFactTol, aDelta);
210     isOnOrIn = Standard_False;
211   }
212
213   aDelta = theU2 - aU2bLast;
214   if(aDelta > aTol)
215   {
216     theU2 = aU2bLast;
217     theFactTol = Max(theFactTol, aDelta);
218     isOnOrIn = Standard_False;
219   }
220
221   aDelta = aV2bFirst - theV2;
222   if(aDelta > aTol)
223   {
224     theV2 = aV2bFirst;
225     theFactTol = Max(theFactTol, aDelta);
226     isOnOrIn = Standard_False;
227   }
228
229   aDelta = theV2 - aV2bLast;
230   if(aDelta > aTol)
231   {
232     theV2 = aV2bLast;
233     theFactTol = Max(theFactTol, aDelta);
234     isOnOrIn = Standard_False;
235   }
236
237   return isOnOrIn;
238 }
239
240 //==================================================================================
241 // function : IntWalk_PWalking::IntWalk_PWalking
242 // purpose  : 
243 //==================================================================================
244 IntWalk_PWalking::IntWalk_PWalking(const Handle(Adaptor3d_HSurface)& Caro1,
245                                    const Handle(Adaptor3d_HSurface)& Caro2,
246                                    const Standard_Real TolTangency,
247                                    const Standard_Real Epsilon,
248                                    const Standard_Real Deflection,
249                                    const Standard_Real Increment ) 
250                                    :
251
252 done(Standard_True),
253 close(Standard_False),
254 fleche(Deflection),
255 tolconf(Epsilon),
256 sensCheminement(1),
257 myIntersectionOn2S(Caro1,Caro2,TolTangency),
258 STATIC_BLOCAGE_SUR_PAS_TROP_GRAND(0),
259 STATIC_PRECEDENT_INFLEXION(0)
260 {
261   Standard_Real KELARG=20.;
262   //
263   pasMax=Increment*0.2; //-- June 25 99 after problems with precision 
264   Um1 = Adaptor3d_HSurfaceTool::FirstUParameter(Caro1);
265   Vm1 = Adaptor3d_HSurfaceTool::FirstVParameter(Caro1);
266   UM1 = Adaptor3d_HSurfaceTool::LastUParameter(Caro1);
267   VM1 = Adaptor3d_HSurfaceTool::LastVParameter(Caro1);
268
269   Um2 = Adaptor3d_HSurfaceTool::FirstUParameter(Caro2);
270   Vm2 = Adaptor3d_HSurfaceTool::FirstVParameter(Caro2);
271   UM2 = Adaptor3d_HSurfaceTool::LastUParameter(Caro2);
272   VM2 = Adaptor3d_HSurfaceTool::LastVParameter(Caro2);
273
274   ResoU1 = Adaptor3d_HSurfaceTool::UResolution(Caro1,Precision::Confusion());
275   ResoV1 = Adaptor3d_HSurfaceTool::VResolution(Caro1,Precision::Confusion());
276
277   ResoU2 = Adaptor3d_HSurfaceTool::UResolution(Caro2,Precision::Confusion());
278   ResoV2 = Adaptor3d_HSurfaceTool::VResolution(Caro2,Precision::Confusion());
279
280   Standard_Real NEWRESO;
281   Standard_Real MAXVAL;
282   Standard_Real MAXVAL2;
283   //
284   MAXVAL  = Abs(Um1);  MAXVAL2 = Abs(UM1);
285   if(MAXVAL2 > MAXVAL) MAXVAL = MAXVAL2;
286   NEWRESO = ResoU1 * MAXVAL ;
287   if(NEWRESO > ResoU1 &&NEWRESO<10) {    ResoU1 = NEWRESO;  }
288
289
290   MAXVAL  = Abs(Um2);   MAXVAL2 = Abs(UM2);
291   if(MAXVAL2 > MAXVAL) MAXVAL = MAXVAL2;
292   NEWRESO = ResoU2 * MAXVAL ;
293   if(NEWRESO > ResoU2 && NEWRESO<10) {     ResoU2 = NEWRESO;  }
294
295
296   MAXVAL  = Abs(Vm1);  MAXVAL2 = Abs(VM1);
297   if(MAXVAL2 > MAXVAL) MAXVAL = MAXVAL2;
298   NEWRESO = ResoV1 * MAXVAL ;
299   if(NEWRESO > ResoV1 && NEWRESO<10) {     ResoV1 = NEWRESO;  }
300
301
302   MAXVAL  = Abs(Vm2);  MAXVAL2 = Abs(VM2);
303   if(MAXVAL2 > MAXVAL) MAXVAL = MAXVAL2;
304   NEWRESO = ResoV2 * MAXVAL ;
305   if(NEWRESO > ResoV2 && NEWRESO<10) {     ResoV2 = NEWRESO;  }
306
307   pasuv[0]=pasMax*Abs(UM1-Um1);
308   pasuv[1]=pasMax*Abs(VM1-Vm1);
309   pasuv[2]=pasMax*Abs(UM2-Um2);
310   pasuv[3]=pasMax*Abs(VM2-Vm2);
311
312   if(ResoU1>0.0001*pasuv[0]) ResoU1=0.00001*pasuv[0];
313   if(ResoV1>0.0001*pasuv[1]) ResoV1=0.00001*pasuv[1];
314   if(ResoU2>0.0001*pasuv[2]) ResoU2=0.00001*pasuv[2];
315   if(ResoV2>0.0001*pasuv[3]) ResoV2=0.00001*pasuv[3];
316
317
318   if(Adaptor3d_HSurfaceTool::IsUPeriodic(Caro1)==Standard_False) { 
319     //UM1+=KELARG*pasuv[0];  Um1-=KELARG*pasuv[0];
320   }
321   else { 
322     Standard_Real t = UM1-Um1; 
323     if(t<Adaptor3d_HSurfaceTool::UPeriod(Caro1)) { 
324       t=0.5*(Adaptor3d_HSurfaceTool::UPeriod(Caro1)-t);
325       t=(t>KELARG*pasuv[0])? KELARG*pasuv[0] : t;
326       UM1+=t;  Um1-=t;
327     }
328   }
329
330   if(Adaptor3d_HSurfaceTool::IsVPeriodic(Caro1)==Standard_False) { 
331     //VM1+=KELARG*pasuv[1];  Vm1-=KELARG*pasuv[1];
332   }
333   else { 
334     Standard_Real t = VM1-Vm1; 
335     if(t<Adaptor3d_HSurfaceTool::VPeriod(Caro1)) { 
336       t=0.5*(Adaptor3d_HSurfaceTool::VPeriod(Caro1)-t);
337       t=(t>KELARG*pasuv[1])? KELARG*pasuv[1] : t;
338       VM1+=t;  Vm1-=t;
339     }
340   }
341
342   if(Adaptor3d_HSurfaceTool::IsUPeriodic(Caro2)==Standard_False) { 
343     //UM2+=KELARG*pasuv[2];  Um2-=KELARG*pasuv[2];
344   }
345   else { 
346     Standard_Real t = UM2-Um2; 
347     if(t<Adaptor3d_HSurfaceTool::UPeriod(Caro2)) { 
348       t=0.5*(Adaptor3d_HSurfaceTool::UPeriod(Caro2)-t);
349       t=(t>KELARG*pasuv[2])? KELARG*pasuv[2] : t;
350       UM2+=t;  Um2-=t;
351     }
352   }
353
354   if(Adaptor3d_HSurfaceTool::IsVPeriodic(Caro2)==Standard_False) {   
355     //VM2+=KELARG*pasuv[3];  Vm2-=KELARG*pasuv[3];
356   }
357   else { 
358     Standard_Real t = VM2-Vm2; 
359     if(t<Adaptor3d_HSurfaceTool::VPeriod(Caro2)) { 
360       t=0.5*(Adaptor3d_HSurfaceTool::VPeriod(Caro2)-t);
361       t=(t>KELARG*pasuv[3])? KELARG*pasuv[3] : t;
362       VM2+=t;  Vm2-=t;
363     }
364   }
365
366   //-- ComputePasInit(pasuv,Um1,UM1,Vm1,VM1,Um2,UM2,Vm2,VM2,Caro1,Caro2);
367
368   for (Standard_Integer i = 0; i<=3;i++) {
369     if(pasuv[i]>10) 
370       pasuv[i] = 10; 
371     pasInit[i] = pasSav[i] = pasuv[i]; 
372   }
373
374
375 }
376 //==================================================================================
377 // function : IntWalk_PWalking
378 // purpose  : 
379 //==================================================================================
380 IntWalk_PWalking::IntWalk_PWalking(const Handle(Adaptor3d_HSurface)& Caro1,
381                                    const Handle(Adaptor3d_HSurface)& Caro2,
382                                    const Standard_Real TolTangency,
383                                    const Standard_Real Epsilon,
384                                    const Standard_Real Deflection,
385                                    const Standard_Real Increment, 
386                                    const Standard_Real U1,
387                                    const Standard_Real V1,
388                                    const Standard_Real U2, 
389                                    const Standard_Real V2)
390                                    :
391
392 done(Standard_True),
393 close(Standard_False),
394 fleche(Deflection),
395 tolconf(Epsilon),
396 sensCheminement(1),
397 myIntersectionOn2S(Caro1,Caro2,TolTangency),
398 STATIC_BLOCAGE_SUR_PAS_TROP_GRAND(0),
399 STATIC_PRECEDENT_INFLEXION(0)
400 {
401   Standard_Real KELARG=20.;
402   //
403   pasMax=Increment*0.2; //-- June 25 99 after problems with precision 
404   //
405   Um1 = Adaptor3d_HSurfaceTool::FirstUParameter(Caro1);
406   Vm1 = Adaptor3d_HSurfaceTool::FirstVParameter(Caro1);
407   UM1 = Adaptor3d_HSurfaceTool::LastUParameter(Caro1);
408   VM1 = Adaptor3d_HSurfaceTool::LastVParameter(Caro1);
409
410   Um2 = Adaptor3d_HSurfaceTool::FirstUParameter(Caro2);
411   Vm2 = Adaptor3d_HSurfaceTool::FirstVParameter(Caro2);
412   UM2 = Adaptor3d_HSurfaceTool::LastUParameter(Caro2);
413   VM2 = Adaptor3d_HSurfaceTool::LastVParameter(Caro2);
414
415   ResoU1 = Adaptor3d_HSurfaceTool::UResolution(Caro1,Precision::Confusion());
416   ResoV1 = Adaptor3d_HSurfaceTool::VResolution(Caro1,Precision::Confusion());
417
418   ResoU2 = Adaptor3d_HSurfaceTool::UResolution(Caro2,Precision::Confusion());
419   ResoV2 = Adaptor3d_HSurfaceTool::VResolution(Caro2,Precision::Confusion());
420   //
421   Standard_Real NEWRESO, MAXVAL, MAXVAL2;
422   //
423   MAXVAL  = Abs(Um1);  
424   MAXVAL2 = Abs(UM1);
425   if(MAXVAL2 > MAXVAL) {
426     MAXVAL = MAXVAL2;
427   }
428   NEWRESO = ResoU1 * MAXVAL ;
429   if(NEWRESO > ResoU1) {
430     ResoU1 = NEWRESO;  
431   }
432   //
433   MAXVAL  = Abs(Um2);   
434   MAXVAL2 = Abs(UM2);
435   if(MAXVAL2 > MAXVAL){
436     MAXVAL = MAXVAL2;
437   }  
438   NEWRESO = ResoU2 * MAXVAL ;
439   if(NEWRESO > ResoU2) {
440     ResoU2 = NEWRESO;  
441   }
442   //
443   MAXVAL  = Abs(Vm1);  
444   MAXVAL2 = Abs(VM1);
445   if(MAXVAL2 > MAXVAL) {
446     MAXVAL = MAXVAL2;
447   }
448   NEWRESO = ResoV1 * MAXVAL ;
449   if(NEWRESO > ResoV1) {    
450     ResoV1 = NEWRESO; 
451   }
452   //
453   MAXVAL  = Abs(Vm2);  
454   MAXVAL2 = Abs(VM2);
455   if(MAXVAL2 > MAXVAL){
456     MAXVAL = MAXVAL2;
457   }  
458   NEWRESO = ResoV2 * MAXVAL ;
459   if(NEWRESO > ResoV2) {  
460     ResoV2 = NEWRESO;
461   }
462   //
463   pasuv[0]=pasMax*Abs(UM1-Um1);
464   pasuv[1]=pasMax*Abs(VM1-Vm1);
465   pasuv[2]=pasMax*Abs(UM2-Um2);
466   pasuv[3]=pasMax*Abs(VM2-Vm2);
467   //
468   if(Adaptor3d_HSurfaceTool::IsUPeriodic(Caro1)==Standard_False) { 
469     UM1+=KELARG*pasuv[0];  
470     Um1-=KELARG*pasuv[0];
471   }
472   else { 
473     Standard_Real t = UM1-Um1; 
474     if(t<Adaptor3d_HSurfaceTool::UPeriod(Caro1)) { 
475       t=0.5*(Adaptor3d_HSurfaceTool::UPeriod(Caro1)-t);
476       t=(t>KELARG*pasuv[0])? KELARG*pasuv[0] : t;
477       UM1+=t;  
478       Um1-=t;
479     }
480   }
481   //
482   if(Adaptor3d_HSurfaceTool::IsVPeriodic(Caro1)==Standard_False) { 
483     VM1+=KELARG*pasuv[1];
484     Vm1-=KELARG*pasuv[1];
485   }
486   else { 
487     Standard_Real t = VM1-Vm1; 
488     if(t<Adaptor3d_HSurfaceTool::VPeriod(Caro1)) { 
489       t=0.5*(Adaptor3d_HSurfaceTool::VPeriod(Caro1)-t);
490       t=(t>KELARG*pasuv[1])? KELARG*pasuv[1] : t;
491       VM1+=t;  Vm1-=t;
492     }
493   }
494   //
495   if(Adaptor3d_HSurfaceTool::IsUPeriodic(Caro2)==Standard_False) { 
496     UM2+=KELARG*pasuv[2];  
497     Um2-=KELARG*pasuv[2];
498   }
499   else { 
500     Standard_Real t = UM2-Um2; 
501     if(t<Adaptor3d_HSurfaceTool::UPeriod(Caro2)) { 
502       t=0.5*(Adaptor3d_HSurfaceTool::UPeriod(Caro2)-t);
503       t=(t>KELARG*pasuv[2])? KELARG*pasuv[2] : t;
504       UM2+=t;  
505       Um2-=t;
506     }
507   }
508
509   if(Adaptor3d_HSurfaceTool::IsVPeriodic(Caro2)==Standard_False) {   
510     VM2+=KELARG*pasuv[3];  
511     Vm2-=KELARG*pasuv[3];
512   }
513   else { 
514     Standard_Real t = VM2-Vm2; 
515     if(t<Adaptor3d_HSurfaceTool::VPeriod(Caro2)) { 
516       t=0.5*(Adaptor3d_HSurfaceTool::VPeriod(Caro2)-t);
517       t=(t>KELARG*pasuv[3])? KELARG*pasuv[3] : t;
518       VM2+=t;  
519       Vm2-=t;
520     }
521   }
522   //-- ComputePasInit(pasuv,Um1,UM1,Vm1,VM1,Um2,UM2,Vm2,VM2,Caro1,Caro2);
523
524   for (Standard_Integer i = 0; i<=3;i++) {
525     pasInit[i] = pasSav[i] = pasuv[i]; 
526   }  
527
528   if(ResoU1>0.0001*pasuv[0]) ResoU1=0.00001*pasuv[0];
529   if(ResoV1>0.0001*pasuv[1]) ResoV1=0.00001*pasuv[1];
530   if(ResoU2>0.0001*pasuv[2]) ResoU2=0.00001*pasuv[2];
531   if(ResoV2>0.0001*pasuv[3]) ResoV2=0.00001*pasuv[3];
532   //
533   TColStd_Array1OfReal Par(1,4);
534   Par(1) = U1;
535   Par(2) = V1;
536   Par(3) = U2;
537   Par(4) = V2;
538   Perform(Par);
539 }
540
541 //==================================================================================
542 // function : PerformFirstPoint
543 // purpose  : 
544 //==================================================================================
545 Standard_Boolean IntWalk_PWalking::PerformFirstPoint  (const TColStd_Array1OfReal& ParDep,
546                                                        IntSurf_PntOn2S& FirstPoint)   
547 {
548   sensCheminement = 1;
549   close = Standard_False;
550   //
551   Standard_Integer i;
552   TColStd_Array1OfReal Param(1,4);
553   //
554   for (i=1; i<=4; ++i) {
555     Param(i) = ParDep(i);
556   }
557   //-- calculate the first solution point
558   math_FunctionSetRoot  Rsnld(myIntersectionOn2S.Function());
559   //
560   myIntersectionOn2S.Perform(Param,Rsnld);
561   if (!myIntersectionOn2S.IsDone())  { 
562     return Standard_False;
563   }
564
565   if (myIntersectionOn2S.IsEmpty()) {
566     return Standard_False;
567   }
568
569   FirstPoint = myIntersectionOn2S.Point();
570   return Standard_True;
571 }
572 //==================================================================================
573 // function : Perform
574 // purpose  : 
575 //==================================================================================
576 void IntWalk_PWalking::Perform(const TColStd_Array1OfReal& ParDep)    
577 {
578   Perform(ParDep,Um1,Vm1,Um2,Vm2,UM1,VM1,UM2,VM2);
579 }
580
581 //=======================================================================
582 //function : SQDistPointSurface
583 //purpose  : Returns square distance between thePnt and theSurf.
584 //            (theU0, theV0) is initial point for extrema
585 //=======================================================================
586 static Standard_Real SQDistPointSurface(const gp_Pnt &thePnt,
587                                         const Adaptor3d_Surface& theSurf,
588                                         const Standard_Real theU0,
589                                         const Standard_Real theV0)
590 {
591   const Extrema_GenLocateExtPS aExtPS(thePnt, theSurf, theU0, theV0,
592                       Precision::PConfusion(), Precision::PConfusion());
593   if(!aExtPS.IsDone())
594     return RealLast();
595   
596   return aExtPS.SquareDistance();
597 }
598
599 //==================================================================================
600 // function : IsTangentExtCheck
601 // purpose  : Additional check if the surfaces are tangent.
602 //            Checks if any point in one surface lie in another surface
603 //            (with given tolerance)
604 //==================================================================================
605 static Standard_Boolean IsTangentExtCheck(const Handle(Adaptor3d_HSurface)& theSurf1,
606                                           const Handle(Adaptor3d_HSurface)& theSurf2,
607                                           const Standard_Real theU10,
608                                           const Standard_Real theV10,
609                                           const Standard_Real theU20,
610                                           const Standard_Real theV20,
611                                           const Standard_Real theArrStep[])
612 {
613   {
614     gp_Pnt aPt;
615     gp_Vec aDu1, aDv1, aDu2, aDv2;
616     theSurf1->D1(theU10, theV10, aPt, aDu1, aDv1);
617     theSurf2->D1(theU20, theV20, aPt, aDu2, aDv2);
618
619     const gp_Vec  aN1(aDu1.Crossed(aDv1)),
620                   aN2(aDu2.Crossed(aDv2));
621     const Standard_Real aDP = aN1.Dot(aN2),
622                         aSQ1 = aN1.SquareMagnitude(),
623                         aSQ2 = aN2.SquareMagnitude();
624
625     if((aSQ1 < RealSmall()) || (aSQ2 < RealSmall()))
626       return Standard_True; //Tangent
627
628     if(aDP*aDP < 0.9998*aSQ1*aSQ2)
629     {//cos(ang N1<->N2) < 0.9999
630       return Standard_False; //Not tangent
631     }
632   }
633
634   const Standard_Real aSQToler = 4.0e-14;
635   const Standard_Integer aNbItems = 4;
636   const Standard_Real aParUS1[aNbItems] = { theU10 + theArrStep[0],
637                                             theU10 - theArrStep[0],
638                                             theU10, theU10};
639   const Standard_Real aParVS1[aNbItems] = { theV10, theV10,
640                                             theV10 + theArrStep[1],
641                                             theV10 - theArrStep[1]};
642   const Standard_Real aParUS2[aNbItems] = { theU20 + theArrStep[2],
643                                             theU20 - theArrStep[2],
644                                             theU20, theU20};
645   const Standard_Real aParVS2[aNbItems] = { theV20, theV20,
646                                             theV20 + theArrStep[3],
647                                             theV20 - theArrStep[3]};
648
649   for(Standard_Integer i = 0; i < aNbItems; i++)
650   {
651     gp_Pnt aP(theSurf1->Value(aParUS1[i], aParVS1[i]));
652     const Standard_Real aSqDist = SQDistPointSurface(aP, theSurf2->Surface(), theU20, theV20);
653     if(aSqDist > aSQToler)
654       return Standard_False;
655   }
656
657   for(Standard_Integer i = 0; i < aNbItems; i++)
658   {
659     gp_Pnt aP(theSurf2->Value(aParUS2[i], aParVS2[i]));
660     const Standard_Real aSqDist = SQDistPointSurface(aP, theSurf1->Surface(), theU10, theV10);
661     if(aSqDist > aSQToler)
662       return Standard_False;
663   }
664
665   return Standard_True;
666 }
667
668 //==================================================================================
669 // function : Perform
670 // purpose  : 
671 //==================================================================================
672 void IntWalk_PWalking::Perform(const TColStd_Array1OfReal& ParDep,
673                                const Standard_Real u1min,
674                                const Standard_Real v1min,
675                                const Standard_Real u2min,
676                                const Standard_Real v2min,
677                                const Standard_Real u1max,
678                                const Standard_Real v1max,
679                                const Standard_Real u2max,
680                                const Standard_Real v2max)
681 {
682   const Standard_Real aSQDistMax = 1.0e-14;
683   //xf
684
685   Standard_Integer NbPasOKConseq=0;
686   TColStd_Array1OfReal Param(1,4);
687   IntImp_ConstIsoparametric ChoixIso;
688   //xt
689   //
690   done = Standard_False;
691   //
692   // Caro1 and Caro2
693   const Handle(Adaptor3d_HSurface)& Caro1 =myIntersectionOn2S.Function().AuxillarSurface1();
694   const Handle(Adaptor3d_HSurface)& Caro2 =myIntersectionOn2S.Function().AuxillarSurface2();
695   //
696   const Standard_Real UFirst1 = Adaptor3d_HSurfaceTool::FirstUParameter(Caro1);
697   const Standard_Real VFirst1 = Adaptor3d_HSurfaceTool::FirstVParameter(Caro1);
698   const Standard_Real ULast1  = Adaptor3d_HSurfaceTool::LastUParameter (Caro1);
699   const Standard_Real VLast1  = Adaptor3d_HSurfaceTool::LastVParameter (Caro1);
700
701   const Standard_Real UFirst2 = Adaptor3d_HSurfaceTool::FirstUParameter(Caro2);
702   const Standard_Real VFirst2 = Adaptor3d_HSurfaceTool::FirstVParameter(Caro2);
703   const Standard_Real ULast2  = Adaptor3d_HSurfaceTool::LastUParameter (Caro2);
704   const Standard_Real VLast2  = Adaptor3d_HSurfaceTool::LastVParameter (Caro2);
705   //
706   ComputePasInit(pasuv,u1min,u1max,v1min,v1max,u2min,u2max,v2min,v2max,
707     Um1,UM1,Vm1,VM1,Um2,UM2,Vm2,VM2,Caro1,Caro2,pasMax+pasMax);
708   //
709   if(pasuv[0]<100.0*ResoU1) {
710     pasuv[0]=100.0*ResoU1; 
711   }
712   if(pasuv[1]<100.0*ResoV1) {
713     pasuv[1]=100.0*ResoV1; 
714   }
715   if(pasuv[2]<100.0*ResoU2) {
716     pasuv[2]=100.0*ResoU2;
717   }
718   if(pasuv[3]<100.0*ResoV2) {
719     pasuv[3]=100.0*ResoV2;
720   }
721   //
722   for (Standard_Integer i=0; i<4; ++i)
723   {
724     if(pasuv[i]>10)
725     {
726       pasuv[i] = 10;
727     }
728
729     pasInit[i] = pasSav[i] = pasuv[i]; 
730   }
731   //
732   line = new IntSurf_LineOn2S ();
733   //
734   for (Standard_Integer i=1; i<=4; ++i)
735   {
736     Param(i)=ParDep(i);
737   }
738   //-- reproduce steps uv connected to surfaces Caro1 and Caro2
739   //-- pasuv[] and pasSav[] are modified during the marching
740   for(Standard_Integer i = 0; i < 4; ++i)
741   {
742     pasSav[i] = pasuv[i] = pasInit[i];
743   }
744
745   //-- calculate the first solution point
746   math_FunctionSetRoot  Rsnld(myIntersectionOn2S.Function());
747   //
748   ChoixIso = myIntersectionOn2S.Perform(Param,Rsnld);
749   if (!myIntersectionOn2S.IsDone())
750   {
751     return;
752   }
753
754   //
755   if (myIntersectionOn2S.IsEmpty())
756   {
757     return;
758   }
759   //
760   if(myIntersectionOn2S.IsTangent())
761   {
762     return;
763   }
764   //
765   Standard_Boolean Arrive, DejaReparti;
766   const Standard_Integer RejectIndexMAX = 250000;
767   Standard_Integer IncKey, RejectIndex;
768   gp_Pnt pf,pl;
769   //
770   DejaReparti = Standard_False;
771   IncKey = 0;
772   RejectIndex = 0;
773   //
774   previousPoint = myIntersectionOn2S.Point();
775   previoustg = Standard_False;
776   previousd  = myIntersectionOn2S.Direction();
777   previousd1 = myIntersectionOn2S.DirectionOnS1();
778   previousd2 = myIntersectionOn2S.DirectionOnS2();
779   indextg = 1;
780   tgdir   = previousd;
781   firstd1 = previousd1;
782   firstd2 = previousd2;
783   tgfirst = tglast = Standard_False;
784   choixIsoSav  =  ChoixIso;
785   //------------------------------------------------------------
786   //-- Test if the first point of marching corresponds 
787   //-- to a point on borders. 
788   //-- In this case, DejaReparti is initialized as True
789   //-- 
790   pf = previousPoint.Value();
791   Standard_Boolean bTestFirstPoint = Standard_True;
792
793   previousPoint.Parameters(Param(1),Param(2),Param(3),Param(4));
794
795   if(IsTangentExtCheck(Caro1, Caro2, Param(1), Param(2), Param(3), Param(4), pasuv))
796     return;
797
798   AddAPoint(line,previousPoint);
799   //
800   IntWalk_StatusDeflection Status = IntWalk_OK;
801   Standard_Boolean NoTestDeflection = Standard_False;
802   Standard_Real SvParam[4], f;
803   Standard_Integer LevelOfEmptyInmyIntersectionOn2S=0;
804   Standard_Integer LevelOfPointConfondu = 0; 
805   Standard_Integer LevelOfIterWithoutAppend = -1;
806   //
807
808   const Standard_Real aTol[4] = { Epsilon(u1max - u1min),
809                                   Epsilon(v1max - v1min),
810                                   Epsilon(u2max - u2min),
811                                   Epsilon(v2max - v2min)};
812   Arrive = Standard_False;
813   while(!Arrive) //010
814   {
815     LevelOfIterWithoutAppend++;
816     if(LevelOfIterWithoutAppend>20)
817     {
818       Arrive = Standard_True; 
819       if(DejaReparti) {
820         break;
821       }
822       RepartirOuDiviser(DejaReparti,ChoixIso,Arrive);
823       LevelOfIterWithoutAppend = 0;
824     }
825     //
826     // compute f
827     f = 0.;
828     switch (ChoixIso) { 
829       case IntImp_UIsoparametricOnCaro1: f = Abs(previousd1.X()); break;
830       case IntImp_VIsoparametricOnCaro1: f = Abs(previousd1.Y()); break;
831       case IntImp_UIsoparametricOnCaro2: f = Abs(previousd2.X()); break;
832       case IntImp_VIsoparametricOnCaro2: f = Abs(previousd2.Y()); break;
833       default:break;
834     }
835     //
836     if(f<0.1) {
837       f=0.1;
838     }
839     //
840     previousPoint.Parameters(Param(1),Param(2),Param(3),Param(4));
841     //
842     //--ofv.begin
843     Standard_Real aIncKey, aEps, dP1, dP2, dP3, dP4;
844     //
845     dP1 = sensCheminement * pasuv[0] * previousd1.X() /f;
846     dP2 = sensCheminement * pasuv[1] * previousd1.Y() /f;
847     dP3 = sensCheminement * pasuv[2] * previousd2.X() /f; 
848     dP4 = sensCheminement * pasuv[3] * previousd2.Y() /f;
849     //
850     aIncKey=5.*(Standard_Real)IncKey;
851     aEps=1.e-7;
852     if(ChoixIso == IntImp_UIsoparametricOnCaro1 && Abs(dP1) < aEps)
853     {
854       dP1 *= aIncKey;
855     }
856
857     if(ChoixIso == IntImp_VIsoparametricOnCaro1 && Abs(dP2) < aEps)
858     {
859       dP2 *= aIncKey;
860     }
861
862     if(ChoixIso == IntImp_UIsoparametricOnCaro2 && Abs(dP3) < aEps)
863     {
864       dP3 *= aIncKey;
865     }
866
867     if(ChoixIso == IntImp_VIsoparametricOnCaro2 && Abs(dP4) < aEps)
868     {
869       dP4 *= aIncKey;
870     }
871     //--ofv.end
872     //
873     Param(1) += dP1;
874     Param(2) += dP2;
875     Param(3) += dP3; 
876     Param(4) += dP4;
877     //==========================
878     SvParam[0]=Param(1); 
879     SvParam[1]=Param(2);
880     SvParam[2]=Param(3);
881     SvParam[3]=Param(4);
882     //
883     Standard_Integer aTryNumber = 0;
884     Standard_Real    isBadPoint = Standard_False;
885     IntImp_ConstIsoparametric aBestIso = ChoixIso;
886     do
887     {
888       isBadPoint = Standard_False;
889
890       ChoixIso= myIntersectionOn2S.Perform(Param, Rsnld, aBestIso);
891
892       if (myIntersectionOn2S.IsDone() && !myIntersectionOn2S.IsEmpty())
893       {
894         //If we go along any surface boundary then it is possible 
895         //to find "outboundaried" point.
896         //Nevertheless, if this deflection is quite small, we will be
897         //able to adjust this point to the boundary.
898
899         Standard_Real aNewPnt[4], anAbsParamDist[4];
900         myIntersectionOn2S.Point().Parameters(aNewPnt[0], aNewPnt[1], aNewPnt[2], aNewPnt[3]);
901         const Standard_Real aParMin[4] = {u1min, v1min, u2min, v2min};
902         const Standard_Real aParMax[4] = {u1max, v1max, u2max, v2max};
903
904         for(Standard_Integer i = 0; i < 4; i++)
905         {
906           if(Abs(aNewPnt[i] - aParMin[i]) < aTol[i])
907             aNewPnt[i] = aParMin[i];
908           else if(Abs(aNewPnt[i] - aParMax[i]) < aTol[i])
909             aNewPnt[i] = aParMax[i];
910         }
911
912         if (aNewPnt[0] < u1min || aNewPnt[0] > u1max ||
913             aNewPnt[1] < v1min || aNewPnt[1] > v1max ||
914             aNewPnt[2] < u2min || aNewPnt[2] > u2max ||
915             aNewPnt[3] < v2min || aNewPnt[3] > v2max)
916         {
917           break; // Out of borders, handle this later.
918         }
919
920         myIntersectionOn2S.ChangePoint().SetValue(aNewPnt[0],
921                                                   aNewPnt[1],
922                                                   aNewPnt[2],
923                                                   aNewPnt[3]);
924
925         anAbsParamDist[0] = Abs(Param(1) - dP1 - aNewPnt[0]);
926         anAbsParamDist[1] = Abs(Param(2) - dP2 - aNewPnt[1]);
927         anAbsParamDist[2] = Abs(Param(3) - dP3 - aNewPnt[2]);
928         anAbsParamDist[3] = Abs(Param(4) - dP4 - aNewPnt[3]);
929         if (anAbsParamDist[0] < ResoU1 &&
930             anAbsParamDist[1] < ResoV1 &&
931             anAbsParamDist[2] < ResoU2 &&
932             anAbsParamDist[3] < ResoV2 &&
933             Status != IntWalk_PasTropGrand)
934         {
935           isBadPoint = Standard_True;
936           aBestIso = IntImp_ConstIsoparametric((aBestIso + 1) % 4);
937         }
938       }
939     } while (isBadPoint && ++aTryNumber <= 4);
940     //
941     if (!myIntersectionOn2S.IsDone())
942     {
943       //end of line, division
944       Arrive = Standard_False;
945       Param(1)=SvParam[0]; 
946       Param(2)=SvParam[1]; 
947       Param(3)=SvParam[2];
948       Param(4)=SvParam[3];
949       RepartirOuDiviser(DejaReparti, ChoixIso, Arrive);
950     }
951     else  //009 
952     {
953       //== Calculation of exact point from Param(.) is possible
954       if (myIntersectionOn2S.IsEmpty())
955       {
956         Standard_Real u1,v1,u2,v2;
957         previousPoint.Parameters(u1,v1,u2,v2);
958         //
959         Arrive = Standard_False;
960         if(u1<UFirst1 || u1>ULast1)
961         {
962           Arrive=Standard_True;
963         }       
964
965         if(u2<UFirst2 || u2>ULast2)
966         {
967           Arrive=Standard_True;
968         }
969
970         if(v1<VFirst1 || v1>VLast1)
971         {
972           Arrive=Standard_True;
973         }
974
975         if(v2<VFirst2 || v2>VLast2)
976         {
977           Arrive=Standard_True;
978         }
979
980         RepartirOuDiviser(DejaReparti,ChoixIso,Arrive);
981         LevelOfEmptyInmyIntersectionOn2S++;
982         //
983         if(LevelOfEmptyInmyIntersectionOn2S>10)
984         {
985           pasuv[0]=pasSav[0]; 
986           pasuv[1]=pasSav[1]; 
987           pasuv[2]=pasSav[2]; 
988           pasuv[3]=pasSav[3];
989         }
990       }
991       else //008
992       {
993         //============================================================
994         //== A point has been found :  T E S T   D E F L E C T I O N 
995         //============================================================
996         if(NoTestDeflection)
997         {
998           NoTestDeflection = Standard_False;
999         }                 
1000         else
1001         {
1002           if(--LevelOfEmptyInmyIntersectionOn2S<=0)
1003           {
1004             LevelOfEmptyInmyIntersectionOn2S=0;
1005             if(LevelOfIterWithoutAppend < 10)
1006             {
1007               Status = TestDeflection();
1008             }                   
1009             else
1010             {
1011               pasuv[0]*=0.5; 
1012               pasuv[1]*=0.5; 
1013               pasuv[2]*=0.5; 
1014               pasuv[3]*=0.5;
1015             }
1016           }
1017         }
1018
1019         //============================================================
1020         //==       T r a i t e m e n t   s u r   S t a t u s        ==
1021         //============================================================
1022         if(LevelOfPointConfondu > 5)
1023         { 
1024           Status = IntWalk_ArretSurPoint; 
1025           LevelOfPointConfondu = 0;  
1026         }
1027         //
1028         if(Status==IntWalk_OK)
1029         { 
1030           NbPasOKConseq++;
1031           if(NbPasOKConseq >= 5)
1032           {
1033             NbPasOKConseq=0;
1034             Standard_Boolean pastroppetit;
1035             Standard_Real t;
1036             //
1037             do
1038             {
1039               pastroppetit=Standard_True;
1040               //
1041               if(pasuv[0]<pasInit[0])
1042               {
1043                 t = (pasInit[0]-pasuv[0])*0.25;
1044                 if(t>0.1*pasInit[0])
1045                 {
1046                   t=0.1*pasuv[0];
1047                 }
1048
1049                 pasuv[0]+=t; 
1050                 pastroppetit=Standard_False;
1051               }
1052
1053               if(pasuv[1]<pasInit[1])
1054               {
1055                 t = (pasInit[1]-pasuv[1])*0.25;
1056                 if(t>0.1*pasInit[1]) {
1057                   t=0.1*pasuv[1];
1058                 }               
1059
1060                 pasuv[1]+=t; 
1061                 pastroppetit=Standard_False;
1062               }
1063
1064               if(pasuv[2]<pasInit[2])
1065               {
1066                 t = (pasInit[2]-pasuv[2])*0.25;
1067                 if(t>0.1*pasInit[2])
1068                 {
1069                   t=0.1*pasuv[2];
1070                 }
1071
1072                 pasuv[2]+=t; 
1073                 pastroppetit=Standard_False;
1074               }
1075
1076               if(pasuv[3]<pasInit[3])
1077               {
1078                 t = (pasInit[3]-pasuv[3])*0.25;
1079                 if(t>0.1*pasInit[3]) {
1080                   t=0.1*pasuv[3];
1081                 }
1082                 pasuv[3]+=t; 
1083                 pastroppetit=Standard_False;
1084               }
1085               if(pastroppetit)
1086               {
1087                 if(pasMax<0.1)
1088                 {
1089                   pasMax*=1.1;
1090                   pasInit[0]*=1.1; 
1091                   pasInit[1]*=1.1; 
1092                   pasInit[2]*=1.1; 
1093                   pasInit[3]*=1.1; 
1094                 }
1095                 else
1096                 {
1097                   pastroppetit=Standard_False;
1098                 }
1099               }
1100             }
1101             while(pastroppetit);
1102           }
1103         }//Status==IntWalk_OK
1104         else
1105           NbPasOKConseq=0;
1106
1107         //
1108         switch(Status)//007 
1109         {
1110         case IntWalk_ArretSurPointPrecedent:
1111           {
1112             Arrive = Standard_False;
1113             RepartirOuDiviser(DejaReparti, ChoixIso, Arrive);
1114             break;
1115           }
1116         case IntWalk_PasTropGrand:
1117           {
1118             Param(1)=SvParam[0];
1119             Param(2)=SvParam[1]; 
1120             Param(3)=SvParam[2]; 
1121             Param(4)=SvParam[3];
1122
1123             if(LevelOfIterWithoutAppend > 5)
1124             {
1125               for (Standard_Integer i = 0; i < 4; i++)
1126               {
1127                 if (pasSav[i] > pasInit[i])
1128                   continue;
1129
1130                 const Standard_Real aDelta = (pasInit[i]-pasSav[i])*0.25;
1131
1132                 if(aDelta > Epsilon(pasInit[i]))
1133                 {
1134                   pasInit[i] -= aDelta;
1135                   LevelOfIterWithoutAppend=0;
1136                 }
1137               }
1138             }
1139
1140             break;
1141           }
1142         case IntWalk_PointConfondu:
1143           {
1144             LevelOfPointConfondu++;
1145
1146             if(LevelOfPointConfondu>5)
1147             {
1148               Standard_Boolean pastroppetit;
1149               //
1150               do
1151               {
1152                 pastroppetit=Standard_True;
1153
1154                 if(pasuv[0]<pasInit[0])
1155                 {
1156                   pasuv[0]+=(pasInit[0]-pasuv[0])*0.25;
1157                   pastroppetit=Standard_False;
1158                 }
1159
1160                 if(pasuv[1]<pasInit[1])
1161                 {
1162                   pasuv[1]+=(pasInit[1]-pasuv[1])*0.25;
1163                   pastroppetit=Standard_False;
1164                 }
1165
1166                 if(pasuv[2]<pasInit[2])
1167                 {
1168                   pasuv[2]+=(pasInit[2]-pasuv[2])*0.25;
1169                   pastroppetit=Standard_False; 
1170                 }
1171
1172                 if(pasuv[3]<pasInit[3])
1173                 {
1174                   pasuv[3]+=(pasInit[3]-pasuv[3])*0.25;
1175                   pastroppetit=Standard_False;
1176                 }
1177
1178                 if(pastroppetit)
1179                 {
1180                   if(pasMax<0.1)
1181                   {
1182                     pasMax*=1.1;
1183                     pasInit[0]*=1.1;
1184                     pasInit[1]*=1.1;
1185                     pasInit[2]*=1.1;
1186                     pasInit[3]*=1.1; 
1187                   }
1188                   else
1189                   {
1190                     pastroppetit=Standard_False;
1191                   }
1192                 }
1193               }
1194               while(pastroppetit);
1195             }
1196
1197             break;
1198           }
1199         case IntWalk_OK:
1200         case IntWalk_ArretSurPoint://006
1201           {
1202             //=======================================================
1203             //== Stop Test t   :  Frame on Param(.)     ==
1204             //=======================================================
1205             //xft arrive here
1206             Arrive = TestArret(DejaReparti,Param,ChoixIso); 
1207             // JMB 30th December 1999. 
1208             // Some statement below should not be put in comment because they are useful.
1209             // See grid CTO 909 A1 which infinitely loops 
1210             if(Arrive==Standard_False && Status==IntWalk_ArretSurPoint)
1211             {
1212               Arrive=Standard_True;
1213 #ifdef OCCT_DEBUG
1214               cout << "IntWalk_PWalking_1.gxx: Problems with intersection"<<endl;
1215 #endif
1216             }
1217
1218             if(Arrive)
1219             {
1220               NbPasOKConseq = -10;
1221             }
1222
1223             if(!Arrive)//005
1224             {
1225               //=====================================================
1226               //== Param(.) is in the limits                       ==
1227               //==  and does not end a closed  line                ==
1228               //=====================================================
1229               //== Check on the current point of myInters
1230               Standard_Boolean pointisvalid = Standard_False;
1231               {
1232                 Standard_Real u1,v1,u2,v2; 
1233                 myIntersectionOn2S.Point().Parameters(u1,v1,u2,v2);
1234
1235                 //
1236                 if(u1 <= UM1  && u2 <= UM2 && v1 <= VM1 && 
1237                   v2 <= VM2  && u1 >= Um1 && u2 >= Um2 &&
1238                   v1 >= Vm1  && v2 >= Vm2)
1239                 {
1240                   pointisvalid=Standard_True;
1241                 }
1242               }
1243
1244               //
1245               if(pointisvalid)
1246               {
1247                 previousPoint = myIntersectionOn2S.Point();
1248                 previoustg = myIntersectionOn2S.IsTangent();
1249
1250                 if(!previoustg)
1251                 {
1252                   previousd  = myIntersectionOn2S.Direction();
1253                   previousd1 = myIntersectionOn2S.DirectionOnS1();
1254                   previousd2 = myIntersectionOn2S.DirectionOnS2();
1255                 }
1256                 //=====================================================
1257                 //== Check on the previous Point
1258                 {
1259                   Standard_Real u1,v1,u2,v2;
1260                   previousPoint.Parameters(u1,v1,u2,v2); 
1261                   if( u1 <= UM1  && u2 <= UM2 && v1 <= VM1 &&
1262                     v2 <= VM2  && u1 >= Um1 && u2 >= Um2 &&
1263                     v1 >= Vm1  && v2 >= Vm2)
1264                   {
1265                     pl = previousPoint.Value();
1266                     if(bTestFirstPoint)
1267                     {
1268                       if(pf.SquareDistance(pl) < aSQDistMax)
1269                       {
1270                         IncKey++;
1271                         if(IncKey == 5000)
1272                           return;
1273                         else
1274                           continue;
1275                       }
1276                       else
1277                       {
1278                         bTestFirstPoint = Standard_False;
1279                       }
1280                     }
1281                     //
1282                     AddAPoint(line,previousPoint);
1283                     RejectIndex++;
1284
1285                     if(RejectIndex >= RejectIndexMAX)
1286                     {
1287                       break;
1288                     }
1289
1290                     //
1291                     LevelOfIterWithoutAppend = 0;
1292                   }
1293                 }
1294               }//pointisvalid
1295               //====================================================
1296
1297               if(Status == IntWalk_ArretSurPoint)
1298               {
1299                 RepartirOuDiviser(DejaReparti,ChoixIso,Arrive);
1300               }
1301               else
1302               {
1303                 if (line->NbPoints() == 2)
1304                 {
1305                   pasSav[0] = pasuv[0];
1306                   pasSav[1] = pasuv[1];
1307                   pasSav[2] = pasuv[2];
1308                   pasSav[3] = pasuv[3];
1309                 }
1310               }
1311             }//005 if(!Arrive)
1312             else  //004
1313             {
1314               if(close)
1315               {
1316                 //================= la ligne est fermee ===============
1317                 AddAPoint(line,line->Value(1)); //ligne fermee
1318                 LevelOfIterWithoutAppend=0;
1319               }
1320               else    //$$$
1321               {
1322                 //====================================================
1323                 //== Param was not in the limits (was reframed)
1324                 //====================================================
1325                 Standard_Boolean bPrevNotTangent = !previoustg || !myIntersectionOn2S.IsTangent();
1326
1327                 IntImp_ConstIsoparametric SauvChoixIso = ChoixIso;
1328                 ChoixIso = myIntersectionOn2S.Perform(Param,Rsnld,ChoixIso);
1329                 //
1330                 if(!myIntersectionOn2S.IsEmpty()) //002
1331                 {
1332                   // mutially outpasses in the square or intersection in corner
1333
1334                   if(TestArret(Standard_True,Param,ChoixIso))
1335                   {
1336                     NbPasOKConseq = -10;
1337                     ChoixIso = myIntersectionOn2S.Perform(Param,Rsnld,ChoixIso);
1338
1339                     if(!myIntersectionOn2S.IsEmpty())
1340                     {
1341                       previousPoint = myIntersectionOn2S.Point();
1342                       previoustg = myIntersectionOn2S.IsTangent();
1343
1344                       if (!previoustg)
1345                       {
1346                         previousd  = myIntersectionOn2S.Direction();
1347                         previousd1 = myIntersectionOn2S.DirectionOnS1();
1348                         previousd2 = myIntersectionOn2S.DirectionOnS2();
1349                       }
1350
1351                       pl = previousPoint.Value();
1352
1353                       if(bTestFirstPoint)
1354                       {
1355                         if(pf.SquareDistance(pl) < aSQDistMax)
1356                         {
1357                           IncKey++;
1358                           if(IncKey == 5000)
1359                             return;
1360                           else
1361                             continue;
1362                         }
1363                         else
1364                         {
1365                           bTestFirstPoint = Standard_False;
1366                         }
1367                       }
1368                       //
1369                       AddAPoint(line,previousPoint);
1370                       RejectIndex++;
1371
1372                       if(RejectIndex >= RejectIndexMAX)
1373                       {
1374                         break;
1375                       }
1376
1377                       //
1378                       LevelOfIterWithoutAppend=0;
1379                       RepartirOuDiviser(DejaReparti,ChoixIso,Arrive);
1380                     }
1381                     else
1382                     {
1383                       //fail framing divides the step
1384                       Arrive = Standard_False;
1385                       RepartirOuDiviser(DejaReparti,ChoixIso,Arrive);
1386                       NoTestDeflection = Standard_True;
1387                       ChoixIso = SauvChoixIso;
1388                     }
1389                   }//if(TestArret())
1390                   else
1391                   {
1392                     // save the last point
1393                     // to revert to it if the current point is out of bounds
1394
1395                     IntSurf_PntOn2S previousPointSave = previousPoint;
1396                     Standard_Boolean previoustgSave   = previoustg;
1397                     gp_Dir previousdSave              = previousd;
1398                     gp_Dir2d previousd1Save           = previousd1;
1399                     gp_Dir2d previousd2Save           = previousd2;
1400
1401                     previousPoint = myIntersectionOn2S.Point();
1402                     previoustg = myIntersectionOn2S.IsTangent();
1403                     Arrive = Standard_False;
1404
1405                     if(!previoustg)
1406                     {
1407                       previousd  = myIntersectionOn2S.Direction();
1408                       previousd1 = myIntersectionOn2S.DirectionOnS1();
1409                       previousd2 = myIntersectionOn2S.DirectionOnS2();
1410                     }
1411
1412                     //========================================
1413                     //== Check on PreviousPoint @@
1414
1415                     {
1416                       Standard_Real u1,v1,u2,v2;
1417                       previousPoint.Parameters(u1,v1,u2,v2);
1418
1419                       //To save initial 2d points
1420                       gp_Pnt2d ParamPntOnS1(Param(1), Param(2));
1421                       gp_Pnt2d ParamPntOnS2(Param(3), Param(4));
1422
1423                       ///////////////////////////
1424                       Param(1) = u1;
1425                       Param(2) = v1;
1426                       Param(3) = u2;
1427                       Param(4) = v2;
1428                       //
1429
1430                       //xf
1431                       Standard_Boolean bFlag1, bFlag2;
1432                       Standard_Real aTol2D=1.e-11;
1433                       //
1434                       bFlag1=u1 >= Um1-aTol2D && v1 >= Vm1-aTol2D && u1 <= UM1+aTol2D && v1 <= VM1+aTol2D;
1435                       bFlag2=u2 >= Um2-aTol2D && v2 >= Vm2-aTol2D && u2 <= UM2+aTol2D && v2 <= VM2+aTol2D;
1436                       if (bFlag1 && bFlag2)
1437                       {
1438                         if (line->NbPoints() > 1)
1439                         {
1440                           IntSurf_PntOn2S prevprevPoint = line->Value(line->NbPoints()-1);
1441                           Standard_Real ppU1, ppV1, ppU2, ppV2;
1442                           prevprevPoint.Parameters(ppU1, ppV1, ppU2, ppV2);
1443                           Standard_Real pU1, pV1, pU2, pV2;
1444                           previousPointSave.Parameters(pU1, pV1, pU2, pV2);
1445                           gp_Vec2d V1onS1(gp_Pnt2d(ppU1, ppV1), gp_Pnt2d(pU1, pV1));
1446                           gp_Vec2d V2onS1(gp_Pnt2d(pU1, pV1), gp_Pnt2d(u1, v1));
1447                           gp_Vec2d V1onS2(gp_Pnt2d(ppU2, ppV2), gp_Pnt2d(pU2, pV2));
1448                           gp_Vec2d V2onS2(gp_Pnt2d(pU2, pV2), gp_Pnt2d(u2, v2));
1449
1450                           const Standard_Real aDot1 = V1onS1 * V2onS1;
1451                           const Standard_Real aDot2 = V1onS2 * V2onS2;
1452
1453                           if ((aDot1 < 0.0) || (aDot2 < 0.0))
1454                           {
1455                             Arrive = Standard_True;
1456                             break;
1457                           }
1458                         }
1459                         /*
1460                         if(u1 <= UM1  && u2 <= UM2 && v1 <= VM1 &&
1461                         v2 <= VM2  && u1 >= Um1 && u2 >= Um2 &&
1462                         v1 >= Vm1  && v2 >= Vm2)  {
1463                         */                      
1464                         //xt
1465                         pl = previousPoint.Value();
1466
1467                         if(bTestFirstPoint)
1468                         {
1469                           if(pf.SquareDistance(pl) < aSQDistMax)
1470                           {
1471                             IncKey++;
1472
1473                             if(IncKey == 5000)
1474                               return;
1475                             else
1476                               continue;
1477                           }
1478                           else
1479                           {
1480                             bTestFirstPoint = Standard_False;
1481                           }
1482                         }
1483
1484                         //To avoid walking around the same point
1485                         //in the tangent zone near a border
1486
1487                         if (previoustg)
1488                         {
1489                           //There are three consecutive points:
1490                           //previousPointSave -> ParamPnt -> curPnt.
1491
1492                           Standard_Real prevU1, prevV1, prevU2, prevV2;
1493                           previousPointSave.Parameters(prevU1, prevV1, prevU2, prevV2);
1494                           gp_Pnt2d prevPntOnS1(prevU1, prevV1), prevPntOnS2(prevU2, prevV2);
1495                           gp_Pnt2d curPntOnS1(u1, v1), curPntOnS2(u2, v2);
1496                           gp_Vec2d PrevToParamOnS1(prevPntOnS1, ParamPntOnS1);
1497                           gp_Vec2d PrevToCurOnS1(prevPntOnS1, curPntOnS1);
1498                           gp_Vec2d PrevToParamOnS2(prevPntOnS2, ParamPntOnS2);
1499                           gp_Vec2d PrevToCurOnS2(prevPntOnS2, curPntOnS2);
1500                           Standard_Real MaxAngle = 3*M_PI/4;
1501                           Standard_Real anAngleS1 = 0.0, anAngleS2 = 0.0;
1502                           const Standard_Real aSQMParS1 = PrevToParamOnS1.SquareMagnitude();
1503                           const Standard_Real aSQMParS2 = PrevToParamOnS2.SquareMagnitude();
1504                           const Standard_Real aSQMCurS1 = PrevToCurOnS1.SquareMagnitude();
1505                           const Standard_Real aSQMCurS2 = PrevToCurOnS2.SquareMagnitude();
1506
1507                           if(aSQMCurS1 < gp::Resolution())
1508                           {
1509                             //We came back to the one of previos point.
1510                             //Therefore, we must break;
1511
1512                             anAngleS1 = M_PI;
1513                           }
1514                           else if(aSQMParS1 < gp::Resolution())
1515                           {
1516                             //We are walking along tangent zone.
1517                             //It should be continued.
1518                             anAngleS1 = 0.0;
1519                           }
1520                           else
1521                           {
1522                             anAngleS1 = Abs(PrevToParamOnS1.Angle(PrevToCurOnS1));
1523                           }
1524
1525                           if(aSQMCurS2 < gp::Resolution())
1526                           {
1527                             //We came back to the one of previos point.
1528                             //Therefore, we must break;
1529
1530                             anAngleS2 = M_PI;
1531                           }
1532                           else if(aSQMParS2 < gp::Resolution())
1533                           {
1534                             //We are walking along tangent zone.
1535                             //It should be continued;
1536                             anAngleS2 = 0.0;
1537                           }
1538                           else
1539                           {
1540                             anAngleS2 = Abs(PrevToParamOnS2.Angle(PrevToCurOnS2));
1541                           }
1542
1543                           if ((anAngleS1 > MaxAngle) && (anAngleS2 > MaxAngle))
1544                           {
1545                             Arrive = Standard_True;
1546                             break;
1547                           }
1548
1549                           {
1550                             //Check singularity.
1551                             //I.e. check if we are walking along direction, which does not
1552                             //result in comming to any point (i.e. derivative 
1553                             //3D-intersection curve along this direction is equal to 0).
1554                             //A sphere with direction {dU=1, dV=0} from point
1555                             //(U=0, V=M_PI/2) can be considered as example for
1556                             //this case (we cannot find another 3D-point if we go thus).
1557
1558                             //Direction chosen along 1st and 2nd surface correspondingly
1559                             const gp_Vec2d  aDirS1(prevPntOnS1, curPntOnS1),
1560                                             aDirS2(prevPntOnS2, curPntOnS2);
1561
1562                             gp_Pnt aPtemp;
1563                             gp_Vec aDuS1, aDvS1, aDuS2, aDvS2;
1564
1565                             myIntersectionOn2S.Function().AuxillarSurface1()->
1566                                   D1(curPntOnS1.X(), curPntOnS1.Y(), aPtemp, aDuS1, aDvS1);
1567                             myIntersectionOn2S.Function().AuxillarSurface2()->
1568                                   D1(curPntOnS2.X(), curPntOnS2.Y(), aPtemp, aDuS2, aDvS2);
1569
1570                             //Derivative WLine along (it is vector-function indeed)
1571                             //directions chosen
1572                             //(https://en.wikipedia.org/wiki/Directional_derivative#Variation_using_only_direction_of_vector).
1573                             //F1 - on the 1st surface, F2 - on the 2nd surface.
1574                             //x, y, z - coordinates of derivative vector.
1575                             const Standard_Real aF1x =  aDuS1.X()*aDirS1.X() + 
1576                                                         aDvS1.X()*aDirS1.Y();
1577                             const Standard_Real aF1y =  aDuS1.Y()*aDirS1.X() +
1578                                                         aDvS1.Y()*aDirS1.Y();
1579                             const Standard_Real aF1z =  aDuS1.Z()*aDirS1.X() +
1580                                                         aDvS1.Z()*aDirS1.Y();
1581                             const Standard_Real aF2x =  aDuS2.X()*aDirS2.X() +
1582                                                         aDvS2.X()*aDirS2.Y();
1583                             const Standard_Real aF2y =  aDuS2.Y()*aDirS2.X() +
1584                                                         aDvS2.Y()*aDirS2.Y();
1585                             const Standard_Real aF2z =  aDuS2.Z()*aDirS2.X() +
1586                                                         aDvS2.Z()*aDirS2.Y();
1587
1588                             const Standard_Real aF1 = aF1x*aF1x + aF1y*aF1y + aF1z*aF1z;
1589                             const Standard_Real aF2 = aF2x*aF2x + aF2y*aF2y + aF2z*aF2z;
1590
1591                             if((aF1 < gp::Resolution()) && (aF2 < gp::Resolution()))
1592                             {
1593                               //All derivative are equal to 0. Therefore, there is
1594                               //no point in going along direction chosen.
1595                               Arrive = Standard_True;
1596                               break;
1597                             }
1598                           }
1599                         }//if (previoustg) cond.
1600
1601                         ////////////////////////////////////////
1602                         AddAPoint(line,previousPoint);
1603                         RejectIndex++;
1604
1605                         if(RejectIndex >= RejectIndexMAX)
1606                         {
1607                           break;
1608                         }
1609
1610                         //
1611
1612                         LevelOfIterWithoutAppend=0;
1613                         Arrive = Standard_True;
1614                       }
1615                       else
1616                       {
1617                         // revert to the last correctly calculated point
1618                         previousPoint = previousPointSave;
1619                         previoustg    = previoustgSave;
1620                         previousd     = previousdSave;
1621                         previousd1    = previousd1Save;
1622                         previousd2    = previousd2Save;
1623                       }
1624                     }
1625
1626                     //
1627                     Standard_Boolean wasExtended = Standard_False;
1628
1629                     if(Arrive && myIntersectionOn2S.IsTangent() && bPrevNotTangent)
1630                     {
1631                       if(ExtendLineInCommonZone(SauvChoixIso, DejaReparti))
1632                       {
1633                         wasExtended = Standard_True;
1634                         Arrive = Standard_False;
1635                         ChoixIso = SauvChoixIso;
1636                       }
1637                     }
1638
1639                     RepartirOuDiviser(DejaReparti,ChoixIso,Arrive);
1640
1641                     if(Arrive && 
1642                       myIntersectionOn2S.IsDone() && !myIntersectionOn2S.IsEmpty() &&
1643                       myIntersectionOn2S.IsTangent() && bPrevNotTangent &&
1644                       !wasExtended)
1645                     {
1646                       if(ExtendLineInCommonZone(SauvChoixIso, DejaReparti))
1647                       {
1648                         wasExtended = Standard_True;
1649                         Arrive = Standard_False;
1650                         ChoixIso = SauvChoixIso;
1651                       }
1652                     }
1653                   }//else !TestArret() $
1654                 }//$$ end successful framing on border (!myIntersectionOn2S.IsEmpty())
1655                 else
1656                 {
1657                   //echec framing on border; division of step 
1658                   Arrive = Standard_False;
1659                   NoTestDeflection = Standard_True;
1660                   RepartirOuDiviser(DejaReparti,ChoixIso,Arrive);
1661                 }
1662               }//$$$ end framing on border (!close)
1663             }//004 fin TestArret return Arrive = True
1664           } // 006case IntWalk_ArretSurPoint:  end Processing Status = OK  or ArretSurPoint 
1665         } //007  switch(Status) 
1666       } //008 end processing point  (TEST DEFLECTION)
1667     } //009 end processing line (else if myIntersectionOn2S.IsDone())
1668   }  //010 end if first departure point allows marching  while (!Arrive)
1669
1670   done = Standard_True;
1671 }
1672 // ===========================================================================================================
1673 // function: ExtendLineInCommonZone
1674 // purpose:  Extends already computed line inside tangent zone in the direction given by theChoixIso.
1675 //           Returns Standard_True if the line was extended through tangent zone and the last computed point 
1676 //           is outside the tangent zone (but it is not put into the line). Otherwise returns Standard_False.
1677 // ===========================================================================================================
1678 Standard_Boolean IntWalk_PWalking::ExtendLineInCommonZone(const IntImp_ConstIsoparametric theChoixIso,
1679                                                           const Standard_Boolean          theDirectionFlag) 
1680 {
1681   Standard_Boolean bOutOfTangentZone = Standard_False;
1682   Standard_Boolean bStop = !myIntersectionOn2S.IsTangent();
1683   Standard_Integer dIncKey = 1;
1684   TColStd_Array1OfReal Param(1,4);
1685   IntWalk_StatusDeflection Status = IntWalk_OK;
1686   Standard_Integer nbIterWithoutAppend = 0;
1687   Standard_Integer nbEqualPoints = 0;
1688   Standard_Integer parit = 0;
1689   Standard_Integer uvit = 0;
1690   IntSurf_SequenceOfPntOn2S aSeqOfNewPoint;
1691
1692   while (!bStop) {
1693     nbIterWithoutAppend++;
1694
1695     if((nbIterWithoutAppend > 20) || (nbEqualPoints > 20)) {
1696 #ifdef OCCT_DEBUG
1697       cout<<"Infinite loop detected. Stop iterations (IntWalk_PWalking_1.gxx)" << endl;
1698 #endif
1699       bStop = Standard_True;
1700       break;
1701     }
1702     Standard_Real f = 0.;
1703
1704     switch (theChoixIso)
1705     { 
1706     case IntImp_UIsoparametricOnCaro1: f = Abs(previousd1.X()); break;
1707     case IntImp_VIsoparametricOnCaro1: f = Abs(previousd1.Y()); break;
1708     case IntImp_UIsoparametricOnCaro2: f = Abs(previousd2.X()); break;
1709     case IntImp_VIsoparametricOnCaro2: f = Abs(previousd2.Y()); break;
1710     }
1711
1712     if(f<0.1) f=0.1;
1713
1714     previousPoint.Parameters(Param(1),Param(2),Param(3),Param(4));
1715
1716     Standard_Real dP1 = sensCheminement * pasuv[0] * previousd1.X() /f;
1717     Standard_Real dP2 = sensCheminement * pasuv[1] * previousd1.Y() /f;
1718     Standard_Real dP3 = sensCheminement * pasuv[2] * previousd2.X() /f; 
1719     Standard_Real dP4 = sensCheminement * pasuv[3] * previousd2.Y() /f;
1720
1721     if(theChoixIso == IntImp_UIsoparametricOnCaro1 && Abs(dP1) < 1.e-7) dP1 *= (5. * (Standard_Real)dIncKey);
1722     if(theChoixIso == IntImp_VIsoparametricOnCaro1 && Abs(dP2) < 1.e-7) dP2 *= (5. * (Standard_Real)dIncKey);
1723     if(theChoixIso == IntImp_UIsoparametricOnCaro2 && Abs(dP3) < 1.e-7) dP3 *= (5. * (Standard_Real)dIncKey);
1724     if(theChoixIso == IntImp_VIsoparametricOnCaro2 && Abs(dP4) < 1.e-7) dP4 *= (5. * (Standard_Real)dIncKey);
1725
1726     Param(1) += dP1;
1727     Param(2) += dP2;
1728     Param(3) += dP3; 
1729     Param(4) += dP4;
1730     Standard_Real SvParam[4];
1731     IntImp_ConstIsoparametric ChoixIso = theChoixIso;
1732
1733     for(parit = 0; parit < 4; parit++) {
1734       SvParam[parit] = Param(parit+1);
1735     }
1736     math_FunctionSetRoot  Rsnld(myIntersectionOn2S.Function());
1737     ChoixIso = myIntersectionOn2S.Perform(Param,Rsnld, theChoixIso);
1738
1739     if (!myIntersectionOn2S.IsDone()) {
1740       return bOutOfTangentZone;
1741     }
1742     else {
1743       if (myIntersectionOn2S.IsEmpty()) {
1744         return bOutOfTangentZone;
1745       }
1746
1747       Status = TestDeflection();
1748
1749       if(Status == IntWalk_OK) {
1750
1751         for(uvit = 0; uvit < 4; uvit++) {
1752           if(pasuv[uvit] < pasInit[uvit]) {
1753             pasuv[uvit] = pasInit[uvit];
1754           }
1755         }
1756       }
1757
1758       switch(Status) {
1759       case  IntWalk_ArretSurPointPrecedent:
1760         {
1761           bStop = Standard_True;
1762           bOutOfTangentZone = !myIntersectionOn2S.IsTangent();
1763           break;
1764         }
1765       case IntWalk_PasTropGrand:
1766         {
1767           for(parit = 0; parit < 4; parit++) {
1768             Param(parit+1) = SvParam[parit];
1769           }
1770           Standard_Boolean bDecrease = Standard_False;
1771
1772           for(uvit = 0; uvit < 4; uvit++) {
1773             if(pasSav[uvit] < pasInit[uvit]) { 
1774               pasInit[uvit] -= (pasInit[uvit] - pasSav[uvit]) * 0.1;
1775               bDecrease = Standard_True;
1776             }
1777           }
1778
1779           if(bDecrease) nbIterWithoutAppend--;
1780           break;
1781         }
1782       case IntWalk_PointConfondu:
1783         {
1784           for(uvit = 0; uvit < 4; uvit++) {
1785             if(pasuv[uvit] < pasInit[uvit]) {
1786               pasuv[uvit] += (pasInit[uvit] - pasuv[uvit]) * 0.1;
1787             }
1788           }
1789           break;
1790         }
1791       case IntWalk_OK:
1792       case IntWalk_ArretSurPoint:
1793         {
1794           //
1795           bStop = TestArret(theDirectionFlag, Param, ChoixIso);
1796           //
1797
1798           //
1799           if(!bStop) {
1800             Standard_Real u11,v11,u12,v12; 
1801             myIntersectionOn2S.Point().Parameters(u11,v11,u12,v12); 
1802             Standard_Real u21,v21,u22,v22;
1803             previousPoint.Parameters(u21,v21,u22,v22); 
1804
1805             if(((fabs(u11-u21) < ResoU1) && (fabs(v11-v21) < ResoV1)) ||
1806               ((fabs(u12-u22) < ResoU2) && (fabs(v12-v22) < ResoV2))) {
1807                 nbEqualPoints++;
1808             }
1809             else {
1810               nbEqualPoints = 0;
1811             }
1812           }
1813           //
1814
1815           bStop = bStop || !myIntersectionOn2S.IsTangent();
1816           bOutOfTangentZone = !myIntersectionOn2S.IsTangent();
1817
1818           if(!bStop) {
1819             Standard_Boolean pointisvalid = Standard_False;
1820             Standard_Real u1,v1,u2,v2; 
1821             myIntersectionOn2S.Point().Parameters(u1,v1,u2,v2); 
1822
1823             if(u1 <= UM1  && u2 <= UM2 && v1 <= VM1 && 
1824               v2 <= VM2  && u1 >= Um1 && u2 >= Um2 &&
1825               v1 >= Vm1  && v2 >= Vm2) 
1826               pointisvalid = Standard_True;
1827
1828             if(pointisvalid) {
1829               previousPoint = myIntersectionOn2S.Point();
1830               previoustg = myIntersectionOn2S.IsTangent();
1831
1832               if(!previoustg) {
1833                 previousd  = myIntersectionOn2S.Direction();
1834                 previousd1 = myIntersectionOn2S.DirectionOnS1();
1835                 previousd2 = myIntersectionOn2S.DirectionOnS2();
1836               }
1837               Standard_Boolean bAddPoint = Standard_True;
1838
1839               if(line->NbPoints() >= 1) {
1840                 gp_Pnt pf = line->Value(1).Value();
1841                 gp_Pnt pl = previousPoint.Value(); 
1842
1843                 if(pf.Distance(pl) < Precision::Confusion()) { 
1844                   dIncKey++; 
1845                   if(dIncKey == 5000) return bOutOfTangentZone; 
1846                   else bAddPoint = Standard_False;
1847                 }
1848               }
1849
1850               if(bAddPoint) {
1851                 aSeqOfNewPoint.Append(previousPoint);
1852                 nbIterWithoutAppend = 0;
1853               }
1854             }
1855
1856             if (line->NbPoints() == 2) {
1857               for(uvit = 0; uvit < 4; uvit++) {
1858                 pasSav[uvit] = pasuv[uvit]; 
1859               }
1860             }
1861
1862             if ( !pointisvalid ) {
1863               // decrease step if out of bounds
1864               // otherwise the same calculations will be 
1865               // repeated several times
1866               if ( ( u1 > UM1 ) || ( u1 < Um1 ) )
1867                 pasuv[0] *= 0.5;
1868
1869               if ( ( v1 > VM1 ) || ( v1 < Vm1 ) ) 
1870                 pasuv[1] *= 0.5;
1871
1872               if ( ( u2 > UM2 ) || ( u2 < Um2 ) )
1873                 pasuv[2] *= 0.5;
1874
1875               if ( ( v2 > VM2 ) || ( v2 < Vm2 ) )
1876                 pasuv[3] *= 0.5;
1877             }
1878           } // end if(!bStop)
1879           else { //if(bStop)
1880             if(close && (line->NbPoints() >= 1)) { 
1881
1882               if(!bOutOfTangentZone) {
1883                 aSeqOfNewPoint.Append(line->Value(1)); // line end
1884               }
1885               nbIterWithoutAppend = 0;
1886             }
1887             else {
1888               ChoixIso = myIntersectionOn2S.Perform(Param, Rsnld, theChoixIso);
1889
1890               if(myIntersectionOn2S.IsEmpty()) { 
1891                 bStop = !myIntersectionOn2S.IsTangent();
1892                 bOutOfTangentZone = !myIntersectionOn2S.IsTangent();
1893               }
1894               else {
1895                 Standard_Boolean bAddPoint = Standard_True;
1896                 Standard_Boolean pointisvalid = Standard_False;
1897
1898                 previousPoint = myIntersectionOn2S.Point();
1899                 Standard_Real u1,v1,u2,v2; 
1900                 previousPoint.Parameters(u1,v1,u2,v2); 
1901
1902                 if(u1 <= UM1  && u2 <= UM2 && v1 <= VM1 && 
1903                   v2 <= VM2  && u1 >= Um1 && u2 >= Um2 &&
1904                   v1 >= Vm1  && v2 >= Vm2) 
1905                   pointisvalid = Standard_True;
1906
1907                 if(pointisvalid) {
1908
1909                   if(line->NbPoints() >= 1) {
1910                     gp_Pnt pf = line->Value(1).Value();
1911                     gp_Pnt pl = previousPoint.Value(); 
1912
1913                     if(pf.Distance(pl) < Precision::Confusion()) { 
1914                       dIncKey++; 
1915                       if(dIncKey == 5000) return bOutOfTangentZone; 
1916                       else bAddPoint = Standard_False;
1917                     }
1918                   }
1919
1920                   if(bAddPoint && !bOutOfTangentZone) {
1921                     aSeqOfNewPoint.Append(previousPoint);
1922                     nbIterWithoutAppend = 0;
1923                   }
1924                 }
1925               }
1926             }
1927           }
1928           break;
1929         }
1930       default:
1931         {
1932           break;
1933         }
1934       }
1935     }
1936   }
1937   Standard_Boolean bExtendLine = Standard_False;
1938   Standard_Real u1 = 0., v1 = 0., u2 = 0., v2 = 0.; 
1939
1940   Standard_Integer pit = 0;
1941
1942   for(pit = 0; !bExtendLine && (pit < 2); pit++) {
1943     if(pit == 0)
1944       previousPoint.Parameters(u1,v1,u2,v2); 
1945     else {
1946       if(aSeqOfNewPoint.Length() > 0)
1947         aSeqOfNewPoint.Value(aSeqOfNewPoint.Length()).Parameters(u1,v1,u2,v2); 
1948       else
1949         break;
1950     }
1951
1952     if(((u1 - Um1) < ResoU1) ||
1953       ((UM1 - u1) < ResoU1) ||
1954       ((u2 - Um2) < ResoU2) ||
1955       ((UM2 - u2) < ResoU2) ||
1956       ((v1 - Vm1) < ResoV1) ||
1957       ((VM1 - v1) < ResoV1) ||
1958       ((v2 - Vm2) < ResoV2) ||
1959       ((VM2 - v2) < ResoV2))
1960       bExtendLine = Standard_True;
1961   }
1962
1963   if(!bExtendLine) {
1964     //    if(Status == IntWalk_OK || Status == IntWalk_ArretSurPoint) {
1965     if(Status == IntWalk_OK) {
1966       bExtendLine = Standard_True;
1967
1968       if(aSeqOfNewPoint.Length() > 1) {
1969         TColStd_Array1OfReal FirstParams(0, 3), LastParams(0, 3), Resolutions(0, 3);
1970         Resolutions(0) = ResoU1; Resolutions(1) = ResoV1; Resolutions(2) = ResoU2; Resolutions(3) = ResoV2;
1971
1972         aSeqOfNewPoint(1).Parameters(FirstParams.ChangeValue(0), FirstParams.ChangeValue(1),
1973           FirstParams.ChangeValue(2), FirstParams.ChangeValue(3));
1974         aSeqOfNewPoint(aSeqOfNewPoint.Length()).Parameters(LastParams.ChangeValue(0), 
1975           LastParams.ChangeValue(1),
1976           LastParams.ChangeValue(2), 
1977           LastParams.ChangeValue(3)); 
1978         Standard_Integer indexofiso = 0;
1979
1980         if(theChoixIso == IntImp_UIsoparametricOnCaro1) indexofiso = 0;
1981         if(theChoixIso == IntImp_VIsoparametricOnCaro1) indexofiso = 1;
1982         if(theChoixIso == IntImp_UIsoparametricOnCaro2) indexofiso = 2;
1983         if(theChoixIso == IntImp_VIsoparametricOnCaro2) indexofiso = 3;
1984
1985         Standard_Integer afirstindex = (indexofiso < 2) ? 0 : 2;
1986         gp_Vec2d aTangentZoneDir(gp_Pnt2d(FirstParams.Value(afirstindex), FirstParams.Value(afirstindex + 1)),
1987           gp_Pnt2d(LastParams.Value(afirstindex), LastParams.Value(afirstindex + 1)));
1988
1989         gp_Dir2d anIsoDir(0, 1);
1990
1991         if((indexofiso == 1) || (indexofiso == 3))
1992           anIsoDir = gp_Dir2d(1, 0);
1993
1994         if(aTangentZoneDir.SquareMagnitude() > gp::Resolution()) {
1995           Standard_Real piquota = M_PI*0.25;
1996
1997           if(fabs(aTangentZoneDir.Angle(anIsoDir)) > piquota) {
1998             Standard_Integer ii = 1, nextii = 2;
1999             gp_Vec2d d1(0, 0);
2000             Standard_Real asqresol = gp::Resolution();
2001             asqresol *= asqresol;
2002
2003             do {
2004               aSeqOfNewPoint(ii).Parameters(FirstParams.ChangeValue(0), FirstParams.ChangeValue(1),
2005                 FirstParams.ChangeValue(2), FirstParams.ChangeValue(3));
2006               aSeqOfNewPoint(ii + 1).Parameters(LastParams.ChangeValue(0), LastParams.ChangeValue(1),
2007                 LastParams.ChangeValue(2), LastParams.ChangeValue(3));
2008               d1 = gp_Vec2d(gp_Pnt2d(FirstParams.Value(afirstindex),
2009                 FirstParams.Value(afirstindex + 1)),
2010                 gp_Pnt2d(LastParams.Value(afirstindex),
2011                 LastParams.Value(afirstindex + 1)));
2012               ii++;
2013             }
2014             while((d1.SquareMagnitude() < asqresol) &&
2015               (ii < aSeqOfNewPoint.Length()));
2016
2017             nextii = ii;
2018
2019             while(nextii < aSeqOfNewPoint.Length()) {
2020
2021               gp_Vec2d nextd1(0, 0);
2022               Standard_Integer jj = nextii;
2023
2024               do {
2025                 aSeqOfNewPoint(jj).Parameters(FirstParams.ChangeValue(0), FirstParams.ChangeValue(1),
2026                   FirstParams.ChangeValue(2), FirstParams.ChangeValue(3));
2027                 aSeqOfNewPoint(jj + 1).Parameters(LastParams.ChangeValue(0), LastParams.ChangeValue(1),
2028                   LastParams.ChangeValue(2), LastParams.ChangeValue(3));
2029                 nextd1 = gp_Vec2d(gp_Pnt2d(FirstParams.Value(afirstindex),
2030                   FirstParams.Value(afirstindex + 1)),
2031                   gp_Pnt2d(LastParams.Value(afirstindex),
2032                   LastParams.Value(afirstindex + 1)));
2033                 jj++;
2034
2035               }
2036               while((nextd1.SquareMagnitude() < asqresol) &&
2037                 (jj < aSeqOfNewPoint.Length()));
2038               nextii = jj;
2039
2040               if(fabs(d1.Angle(nextd1)) > piquota) {
2041                 bExtendLine = Standard_False;
2042                 break;
2043               }
2044               d1 = nextd1;
2045             }
2046           }
2047           // end if(fabs(aTangentZoneDir.Angle(anIsoDir)
2048         }
2049       }
2050     }
2051   }
2052
2053   if(!bExtendLine) {
2054     return Standard_False;
2055   }
2056   Standard_Integer i = 0;
2057
2058   for(i = 1; i <= aSeqOfNewPoint.Length(); i++) {
2059     AddAPoint(line, aSeqOfNewPoint.Value(i));
2060   }
2061
2062   return bOutOfTangentZone;
2063 }
2064
2065 //=======================================================================
2066 //function : DistanceMinimizeByGradient
2067 //purpose  : 
2068 //=======================================================================
2069 Standard_Boolean IntWalk_PWalking::
2070 DistanceMinimizeByGradient( const Handle(Adaptor3d_HSurface)& theASurf1,
2071                            const Handle(Adaptor3d_HSurface)& theASurf2,
2072                            Standard_Real& theU1,
2073                            Standard_Real& theV1,
2074                            Standard_Real& theU2,
2075                            Standard_Real& theV2,
2076                            const Standard_Real theStep0U1V1,
2077                            const Standard_Real theStep0U2V2)
2078 {
2079   const Standard_Integer aNbIterMAX = 60;
2080   const Standard_Real aTol = 1.0e-14;
2081   Handle(Geom_Surface) aS1, aS2;
2082
2083   if (theASurf1->GetType() != GeomAbs_BezierSurface &&
2084       theASurf1->GetType() != GeomAbs_BSplineSurface)
2085       return Standard_True;
2086   if (theASurf2->GetType() != GeomAbs_BezierSurface &&
2087       theASurf2->GetType() != GeomAbs_BSplineSurface)
2088       return Standard_True;
2089
2090   Standard_Boolean aStatus = Standard_False;
2091
2092   gp_Pnt aP1, aP2;
2093   gp_Vec aD1u, aD1v, aD2U, aD2V;
2094
2095   theASurf1->D1(theU1, theV1, aP1, aD1u, aD1v);
2096   theASurf2->D1(theU2, theV2, aP2, aD2U, aD2V);
2097
2098   Standard_Real aSQDistPrev = aP1.SquareDistance(aP2);
2099
2100   gp_Vec aP12(aP1, aP2);
2101
2102   Standard_Real aGradFu(-aP12.Dot(aD1u));
2103   Standard_Real aGradFv(-aP12.Dot(aD1v));
2104   Standard_Real aGradFU( aP12.Dot(aD2U));
2105   Standard_Real aGradFV( aP12.Dot(aD2V));
2106
2107   Standard_Real aSTEPuv = theStep0U1V1, aStepUV = theStep0U2V2;
2108
2109   Standard_Boolean flRepeat = Standard_True;
2110   Standard_Integer aNbIter = aNbIterMAX;
2111
2112   while(flRepeat)
2113   {
2114     Standard_Real anAdd = aGradFu*aSTEPuv;
2115     Standard_Real aPARu = (anAdd >= 0.0)?
2116       (theU1 - Max(anAdd, Epsilon(theU1))) :
2117     (theU1 + Max(-anAdd, Epsilon(theU1)));
2118     anAdd = aGradFv*aSTEPuv;
2119     Standard_Real aPARv = (anAdd >= 0.0)?
2120       (theV1 - Max(anAdd, Epsilon(theV1))) :
2121     (theV1 + Max(-anAdd, Epsilon(theV1)));
2122     anAdd = aGradFU*aStepUV;
2123     Standard_Real aParU = (anAdd >= 0.0)?
2124       (theU2 - Max(anAdd, Epsilon(theU2))) :
2125     (theU2 + Max(-anAdd, Epsilon(theU2)));
2126     anAdd = aGradFV*aStepUV;
2127     Standard_Real aParV = (anAdd >= 0.0)?
2128       (theV2 - Max(anAdd, Epsilon(theV2))) :
2129     (theV2 + Max(-anAdd, Epsilon(theV2)));
2130
2131     gp_Pnt aPt1, aPt2;
2132
2133     theASurf1->D1(aPARu, aPARv, aPt1, aD1u, aD1v);
2134     theASurf2->D1(aParU, aParV, aPt2, aD2U, aD2V);
2135
2136     Standard_Real aSQDist = aPt1.SquareDistance(aPt2);
2137
2138     if(aSQDist < aSQDistPrev)
2139     {
2140       aSQDistPrev = aSQDist;
2141       theU1 = aPARu;
2142       theV1 = aPARv;
2143       theU2 = aParU;
2144       theV2 = aParV;
2145
2146       aStatus = aSQDistPrev < aTol;
2147       aSTEPuv *= 1.2;
2148       aStepUV *= 1.2;
2149     }
2150     else
2151     {
2152       if(--aNbIter < 0)
2153       {
2154         flRepeat = Standard_False;
2155       }
2156       else
2157       {
2158         theASurf1->D1(theU1, theV1, aPt1, aD1u, aD1v);
2159         theASurf2->D1(theU2, theV2, aPt2, aD2U, aD2V);
2160
2161         gp_Vec aP12(aPt1, aPt2);
2162         aGradFu = -aP12.Dot(aD1u);
2163         aGradFv = -aP12.Dot(aD1v);
2164         aGradFU = aP12.Dot(aD2U);
2165         aGradFV = aP12.Dot(aD2V);
2166         aSTEPuv = theStep0U1V1;
2167         aStepUV = theStep0U2V2;
2168       }
2169     }
2170   }
2171
2172   return aStatus;
2173 }
2174
2175 //=======================================================================
2176 //function : DistanceMinimizeByExtrema
2177 //purpose  : 
2178 //=======================================================================
2179 Standard_Boolean IntWalk_PWalking::
2180 DistanceMinimizeByExtrema(const Handle(Adaptor3d_HSurface)& theASurf, 
2181                           const gp_Pnt& theP0,
2182                           Standard_Real& theU0,
2183                           Standard_Real& theV0,
2184                           const Standard_Real theStep0U,
2185                           const Standard_Real theStep0V)
2186 {
2187   const Standard_Real aTol = 1.0e-14;
2188   gp_Pnt aPS;
2189   gp_Vec aD1Su, aD1Sv, aD2Su, aD2Sv, aD2SuvTemp;
2190   Standard_Real aSQDistPrev = RealLast();
2191   Standard_Real aU = theU0, aV = theV0;
2192
2193   Standard_Integer aNbIter = 10;
2194   do
2195   {
2196     theASurf->D2(aU, aV, aPS, aD1Su, aD1Sv, aD2Su, aD2Sv, aD2SuvTemp);
2197
2198     gp_Vec aVec(theP0, aPS);
2199
2200     Standard_Real aSQDist = aVec.SquareMagnitude();
2201
2202     if(aSQDist >= aSQDistPrev)
2203       break;
2204
2205     aSQDistPrev = aSQDist;
2206     theU0 = aU;
2207     theV0 = aV;
2208     aNbIter--;
2209
2210     if(aSQDistPrev < aTol)
2211       break;
2212
2213     //Functions
2214     const Standard_Real aF1 = aD1Su.Dot(aVec), aF2 = aD1Sv.Dot(aVec);
2215
2216     //Derivatives
2217     const Standard_Real aDf1u = aD2Su.Dot(aVec) + aD1Su.Dot(aD1Su),
2218       aDf1v = aD2Su.Dot(aD1Sv),
2219       aDf2u = aDf1v,
2220       aDf2v = aD2Sv.Dot(aVec) + aD1Sv.Dot(aD1Sv);
2221
2222     const Standard_Real aDet = aDf1u*aDf2v - aDf1v*aDf2u;
2223     aU -= theStep0U*(aDf2v*aF1 - aDf1v*aF2)/aDet;
2224     aV += theStep0V*(aDf2u*aF1 - aDf1u*aF2)/aDet;
2225   }
2226   while(aNbIter > 0);
2227
2228   return (aSQDistPrev < aTol);
2229 }
2230
2231 //=======================================================================
2232 //function : SeekPointOnBoundary
2233 //purpose  : 
2234 //=======================================================================
2235 Standard_Boolean IntWalk_PWalking::
2236 SeekPointOnBoundary(const Handle(Adaptor3d_HSurface)& theASurf1,
2237                     const Handle(Adaptor3d_HSurface)& theASurf2,
2238                     const Standard_Real theU1,
2239                     const Standard_Real theV1,
2240                     const Standard_Real theU2,
2241                     const Standard_Real theV2,
2242                     const Standard_Boolean isTheFirst)
2243 {
2244   const Standard_Real aTol = 1.0e-14;
2245   Standard_Boolean isOK = Standard_False;
2246   Standard_Real U1prec = theU1, V1prec = theV1, U2prec = theU2, V2prec = theV2;
2247
2248   Standard_Boolean flFinish = Standard_False;
2249
2250   Standard_Integer aNbIter = 20;
2251   while(!flFinish)
2252   {
2253     flFinish = Standard_False;
2254     Standard_Boolean aStatus = Standard_False;
2255
2256     do
2257     {
2258       aNbIter--;
2259       aStatus = DistanceMinimizeByGradient(theASurf1, theASurf2, U1prec, V1prec, U2prec, V2prec);
2260       if(aStatus)
2261       {
2262         break;
2263       }
2264
2265       aStatus = DistanceMinimizeByExtrema(theASurf1, theASurf2->Value(U2prec, V2prec), U1prec, V1prec);
2266       if(aStatus)
2267       {
2268         break;
2269       }
2270
2271       aStatus = DistanceMinimizeByExtrema(theASurf2, theASurf1->Value(U1prec, V1prec), U2prec, V2prec);
2272       if(aStatus)
2273       {
2274         break;
2275       }
2276     }
2277     while(!aStatus && (aNbIter > 0));
2278
2279     if(aStatus)
2280     {
2281       const Standard_Real aTolMax = 1.0e-8;
2282       Standard_Real aTolF = 0.0;
2283
2284       Standard_Real u1 = U1prec, v1 = V1prec, u2 = U2prec, v2 = V2prec;
2285
2286       flFinish = Checking(theASurf1, theASurf2, U1prec, V1prec, U2prec, V2prec, aTolF);
2287
2288       if(aTolF <= aTolMax)
2289       {
2290         gp_Pnt  aP1 = theASurf1->Value(u1, v1),
2291           aP2 = theASurf2->Value(u2, v2);
2292         gp_Pnt aPInt(0.5*(aP1.XYZ() + aP2.XYZ()));
2293
2294         const Standard_Real aSQDist1 = aPInt.SquareDistance(aP1),
2295           aSQDist2 = aPInt.SquareDistance(aP2);
2296         if((aSQDist1 < aTol) && (aSQDist2 < aTol))
2297         {
2298           IntSurf_PntOn2S anIP;
2299           anIP.SetValue(aPInt, u1, v1, u2, v2);
2300
2301           if(isTheFirst)
2302             line->InsertBefore(1,anIP);
2303           else
2304             line->Add(anIP);
2305
2306           isOK = Standard_True;
2307         }
2308       }
2309     }
2310     else
2311     {
2312       break;
2313     }
2314
2315     if(aNbIter < 0)
2316       break;
2317   }
2318
2319   return isOK;
2320 }
2321
2322 //=======================================================================
2323 //function : PutToBoundary
2324 //purpose  : 
2325 //=======================================================================
2326 Standard_Boolean IntWalk_PWalking::
2327 PutToBoundary(const Handle(Adaptor3d_HSurface)& theASurf1,
2328               const Handle(Adaptor3d_HSurface)& theASurf2)
2329 {
2330   const Standard_Real aTolMin = Precision::Confusion();
2331
2332   Standard_Boolean hasBeenAdded = Standard_False;
2333
2334   const Standard_Real aU1bFirst = theASurf1->FirstUParameter();
2335   const Standard_Real aU1bLast = theASurf1->LastUParameter();
2336   const Standard_Real aU2bFirst = theASurf2->FirstUParameter();
2337   const Standard_Real aU2bLast = theASurf2->LastUParameter();
2338   const Standard_Real aV1bFirst = theASurf1->FirstVParameter();
2339   const Standard_Real aV1bLast = theASurf1->LastVParameter();
2340   const Standard_Real aV2bFirst = theASurf2->FirstVParameter();
2341   const Standard_Real aV2bLast = theASurf2->LastVParameter();
2342
2343   Standard_Real aTol = 1.0;
2344   aTol = Min(aTol, aU1bLast - aU1bFirst);
2345   aTol = Min(aTol, aU2bLast - aU2bFirst);
2346   aTol = Min(aTol, aV1bLast - aV1bFirst);
2347   aTol = Min(aTol, aV2bLast - aV2bFirst)*1.0e-3;
2348
2349   if(aTol <= 2.0*aTolMin)
2350     return hasBeenAdded;
2351
2352   Standard_Boolean isNeedAdding = Standard_False;
2353   Standard_Boolean isU1parallel = Standard_False, isV1parallel = Standard_False;
2354   Standard_Boolean isU2parallel = Standard_False, isV2parallel = Standard_False;
2355   IsParallel(line, Standard_True, aTol, isU1parallel, isV1parallel);
2356   IsParallel(line, Standard_False, aTol, isU2parallel, isV2parallel);
2357
2358   Standard_Real u1, v1, u2, v2;
2359   line->Value(1).Parameters(u1, v1, u2, v2);
2360   Standard_Real aDelta = 0.0;
2361
2362   if(!isV1parallel)
2363   {
2364     aDelta = u1 - aU1bFirst;
2365     if((aTolMin < aDelta) && (aDelta < aTol))
2366     {
2367       u1 = aU1bFirst - aDelta;
2368       isNeedAdding = Standard_True;
2369     }
2370     else
2371     {
2372       aDelta = aU1bLast - u1;
2373       if((aTolMin < aDelta) && (aDelta < aTol))
2374       {
2375         u1 = aU1bLast + aDelta;
2376         isNeedAdding = Standard_True;
2377       }
2378     }
2379   }
2380
2381   if(!isV2parallel)
2382   {
2383     aDelta = u2 - aU2bFirst;
2384     if((aTolMin < aDelta) && (aDelta < aTol))
2385     {
2386       u2 = aU2bFirst - aDelta;
2387       isNeedAdding = Standard_True;
2388     }
2389     else
2390     {
2391       aDelta = aU2bLast - u2;
2392       if((aTolMin < aDelta) && (aDelta < aTol))
2393       {
2394         u2 = aU2bLast + aDelta;
2395         isNeedAdding = Standard_True;
2396       }
2397     }
2398   }
2399
2400   if(!isU1parallel)
2401   {
2402     aDelta = v1 - aV1bFirst;
2403     if((aTolMin < aDelta) && (aDelta < aTol))
2404     {
2405       v1 = aV1bFirst - aDelta;
2406       isNeedAdding = Standard_True;
2407     }
2408     else
2409     {
2410       aDelta = aV1bLast - v1;
2411       if((aTolMin < aDelta) && (aDelta < aTol))
2412       {
2413         v1 = aV1bLast + aDelta;
2414         isNeedAdding = Standard_True;
2415       }
2416     }
2417   }
2418
2419   if(!isU2parallel)
2420   {
2421     aDelta = v2 - aV2bFirst;
2422     if((aTolMin < aDelta) && (aDelta < aTol))
2423     {
2424       v2 = aV2bFirst - aDelta;
2425       isNeedAdding = Standard_True;
2426     }
2427     else
2428     {
2429       aDelta = aV2bLast - v2;
2430       if((aTolMin < aDelta) && (aDelta < aTol))
2431       {
2432         v2 = aV2bLast + aDelta;
2433         isNeedAdding = Standard_True;
2434       }
2435     }
2436   }
2437
2438   if(isNeedAdding)
2439   {
2440     hasBeenAdded = 
2441       SeekPointOnBoundary(theASurf1, theASurf2, u1, 
2442       v1, u2, v2, Standard_True);
2443   }
2444
2445   const Standard_Integer aNbPnts = line->NbPoints();
2446   isNeedAdding = Standard_False;
2447   line->Value(aNbPnts).Parameters(u1, v1, u2, v2);
2448
2449   if(!isV1parallel)
2450   {
2451     aDelta = u1 - aU1bFirst;
2452     if((aTolMin < aDelta) && (aDelta < aTol))
2453     {
2454       u1 = aU1bFirst - aDelta;
2455       isNeedAdding = Standard_True;
2456     }
2457     else
2458     {
2459       aDelta = aU1bLast - u1;
2460       if((aTolMin < aDelta) && (aDelta < aTol))
2461       {
2462         u1 = aU1bLast + aDelta;
2463         isNeedAdding = Standard_True;
2464       }
2465     }
2466   }
2467
2468   if(!isV2parallel)
2469   {
2470     aDelta = u2 - aU2bFirst;
2471     if((aTolMin < aDelta) && (aDelta < aTol))
2472     {
2473       u2 = aU2bFirst - aDelta;
2474       isNeedAdding = Standard_True;
2475     }
2476     else
2477     {
2478       aDelta = aU2bLast - u2;
2479       if((aTolMin < aDelta) && (aDelta < aTol))
2480       {
2481         u2 = aU2bLast + aDelta;
2482         isNeedAdding = Standard_True;
2483       }
2484     }
2485   }
2486
2487   if(!isU1parallel)
2488   {
2489     aDelta = v1 - aV1bFirst;
2490     if((aTolMin < aDelta) && (aDelta < aTol))
2491     {
2492       v1 = aV1bFirst - aDelta;
2493       isNeedAdding = Standard_True;
2494     }
2495     else
2496     {
2497       aDelta = aV1bLast - v1;
2498       if((aTolMin < aDelta) && (aDelta < aTol))
2499       {
2500         v1 = aV1bLast + aDelta;
2501         isNeedAdding = Standard_True;
2502       }
2503     }
2504   }
2505
2506   if(!isU2parallel)
2507   {
2508     aDelta = v2 - aV2bFirst;
2509     if((aTolMin < aDelta) && (aDelta < aTol))
2510     {
2511       v2 = aV2bFirst - aDelta;
2512       isNeedAdding = Standard_True;
2513     }
2514     else
2515     {
2516       aDelta = aV2bLast - v2;
2517       if((aTolMin < aDelta) && (aDelta < aTol))
2518       {
2519         v2 = aV2bLast + aDelta;
2520         isNeedAdding = Standard_True;
2521       }
2522     }
2523   }
2524
2525   if(isNeedAdding)
2526   {
2527     hasBeenAdded = 
2528       SeekPointOnBoundary(theASurf1, theASurf2, u1, 
2529       v1, u2, v2, Standard_False);
2530   }
2531
2532   return hasBeenAdded;
2533 }
2534
2535 //=======================================================================
2536 //function : SeekAdditionalPoints
2537 //purpose  : 
2538 //=======================================================================
2539 Standard_Boolean IntWalk_PWalking::
2540 SeekAdditionalPoints( const Handle(Adaptor3d_HSurface)& theASurf1,
2541                      const Handle(Adaptor3d_HSurface)& theASurf2,
2542                      const Standard_Integer theMinNbPoints)
2543 {
2544   const Standard_Real aTol = 1.0e-14;
2545   Standard_Integer aNbPoints = line->NbPoints();
2546   if(aNbPoints > theMinNbPoints)
2547     return Standard_True;
2548
2549   const Standard_Real aU1bFirst = theASurf1->FirstUParameter();
2550   const Standard_Real aU1bLast = theASurf1->LastUParameter();
2551   const Standard_Real aU2bFirst = theASurf2->FirstUParameter();
2552   const Standard_Real aU2bLast = theASurf2->LastUParameter();
2553   const Standard_Real aV1bFirst = theASurf1->FirstVParameter();
2554   const Standard_Real aV1bLast = theASurf1->LastVParameter();
2555   const Standard_Real aV2bFirst = theASurf2->FirstVParameter();
2556   const Standard_Real aV2bLast = theASurf2->LastVParameter();
2557
2558
2559   Standard_Boolean isPrecise = Standard_False;
2560
2561   Standard_Real U1prec = 0.0, V1prec = 0.0, U2prec = 0.0, V2prec = 0.0;
2562
2563   Standard_Integer aNbPointsPrev = 0;
2564   while(aNbPoints < theMinNbPoints && (aNbPoints != aNbPointsPrev))
2565   {
2566     aNbPointsPrev = aNbPoints;
2567     for(Standard_Integer fp = 1, lp = 2; fp < aNbPoints; fp = lp + 1)
2568     {
2569       Standard_Real U1f, V1f, U2f, V2f; //first point in 1st and 2nd surafaces
2570       Standard_Real U1l, V1l, U2l, V2l; //last  point in 1st and 2nd surafaces
2571
2572       lp = fp+1;
2573       line->Value(fp).Parameters(U1f, V1f, U2f, V2f);
2574       line->Value(lp).Parameters(U1l, V1l, U2l, V2l);
2575
2576       U1prec = 0.5*(U1f+U1l);
2577       if(U1prec < aU1bFirst)
2578         U1prec = aU1bFirst;
2579       if(U1prec > aU1bLast)
2580         U1prec = aU1bLast;
2581
2582       V1prec = 0.5*(V1f+V1l);
2583       if(V1prec < aV1bFirst)
2584         V1prec = aV1bFirst;
2585       if(V1prec > aV1bLast)
2586         V1prec = aV1bLast;
2587
2588       U2prec = 0.5*(U2f+U2l);
2589       if(U2prec < aU2bFirst)
2590         U2prec = aU2bFirst;
2591       if(U2prec > aU2bLast)
2592         U2prec = aU2bLast;
2593
2594       V2prec = 0.5*(V2f+V2l);
2595       if(V2prec < aV2bFirst)
2596         V2prec = aV2bFirst;
2597       if(V2prec > aV2bLast)
2598         V2prec = aV2bLast;
2599
2600       Standard_Boolean aStatus = Standard_False;
2601       Standard_Integer aNbIter = 5;
2602       do
2603       {
2604         aStatus = DistanceMinimizeByGradient(theASurf1, theASurf2, U1prec, V1prec, U2prec, V2prec);
2605         if(aStatus)
2606         {
2607           break;
2608         }
2609
2610         aStatus = DistanceMinimizeByExtrema(theASurf1, theASurf2->Value(U2prec, V2prec), U1prec, V1prec);
2611         if(aStatus)
2612         {
2613           break;
2614         }
2615
2616         aStatus = DistanceMinimizeByExtrema(theASurf2, theASurf1->Value(U1prec, V1prec), U2prec, V2prec);
2617         if(aStatus)
2618         {
2619           break;
2620         }
2621       }
2622       while(!aStatus && (--aNbIter > 0));
2623
2624       if(aStatus)
2625       {
2626         gp_Pnt  aP1 = theASurf1->Value(U1prec, V1prec),
2627           aP2 = theASurf2->Value(U2prec, V2prec);
2628         gp_Pnt aPInt(0.5*(aP1.XYZ() + aP2.XYZ()));
2629
2630         const Standard_Real aSQDist1 = aPInt.SquareDistance(aP1),
2631           aSQDist2 = aPInt.SquareDistance(aP2);
2632
2633         if((aSQDist1 < aTol) && (aSQDist2 < aTol))
2634         {
2635           IntSurf_PntOn2S anIP;
2636           anIP.SetValue(aPInt, U1prec, V1prec, U2prec, V2prec);
2637           line->InsertBefore(lp, anIP);
2638
2639           isPrecise = Standard_True;
2640
2641           if(++aNbPoints >= theMinNbPoints)
2642             break;
2643         }
2644         else
2645         {
2646           lp--;
2647         }
2648       }
2649     }
2650   }
2651
2652   return isPrecise;
2653 }
2654
2655 void IntWalk_PWalking::
2656 RepartirOuDiviser(Standard_Boolean& DejaReparti,
2657                   IntImp_ConstIsoparametric& ChoixIso,
2658                   Standard_Boolean& Arrive) 
2659
2660                   // at the neighborhood of a point, there is a fail of marching 
2661                   // it is required to divide the steps to try to continue
2662                   // if the step is too small if we are on border
2663                   // restart in another direction if it was not done, otherwise stop
2664
2665 {
2666   //  Standard_Integer i;
2667   if (Arrive) {    //restart in the other direction
2668     if (!DejaReparti ) {
2669       Arrive        = Standard_False; 
2670       DejaReparti   = Standard_True;
2671       previousPoint = line->Value(1);
2672       previoustg    = Standard_False;
2673       previousd1    = firstd1;
2674       previousd2    = firstd2;
2675       previousd     = tgdir;
2676       indextg       = line->NbPoints();
2677       tgdir.Reverse();
2678       line->Reverse();
2679
2680       //-- printf("\nIntWalk_PWalking_2.gxx Reverse %3d\n",indextg);
2681       sensCheminement = -1;
2682       tgfirst      = tglast;
2683       tglast       = Standard_False;
2684       ChoixIso     = choixIsoSav;
2685 #if 0
2686       pasuv[0]=pasSav[0];
2687       pasuv[1]=pasSav[1];
2688       pasuv[2]=pasSav[2];
2689       pasuv[3]=pasSav[3];
2690 #else 
2691       Standard_Real u1,v1,u2,v2;
2692       Standard_Real U1,V1,U2,V2;
2693       Standard_Integer nn=line->NbPoints();
2694       if(nn>2) { 
2695         line->Value(nn).Parameters(u1,v1,u2,v2);
2696         line->Value(nn-1).Parameters(U1,V1,U2,V2);
2697         pasuv[0]=Abs(u1-U1);
2698         pasuv[1]=Abs(v1-V1);
2699         pasuv[2]=Abs(u2-U2);
2700         pasuv[3]=Abs(v2-V2);
2701       }
2702 #endif
2703
2704     }
2705   }  
2706   else  {
2707     if (    pasuv[0]*0.5 < ResoU1
2708       &&  pasuv[1]*0.5 < ResoV1
2709       &&  pasuv[2]*0.5 < ResoU2
2710       &&  pasuv[3]*0.5 < ResoV2
2711       ) {
2712         if (!previoustg) {
2713           tglast = Standard_True;      // IS IT ENOUGH ????
2714         }
2715
2716         if (!DejaReparti) {  //restart in the other direction
2717           DejaReparti       = Standard_True;
2718           previousPoint     = line->Value(1);
2719           previoustg        = Standard_False;
2720           previousd1        = firstd1;
2721           previousd2        = firstd2;
2722           previousd         = tgdir;
2723           indextg           = line->NbPoints();
2724           tgdir.Reverse();
2725           line->Reverse();
2726
2727           //-- printf("\nIntWalk_PWalking_2.gxx Reverse %3d\n",indextg);
2728
2729           sensCheminement   = -1;
2730           tgfirst           = tglast;
2731           tglast            = Standard_False;
2732           ChoixIso          = choixIsoSav;
2733
2734 #if 0 
2735           pasuv[0]=pasSav[0];
2736           pasuv[1]=pasSav[1];
2737           pasuv[2]=pasSav[2];
2738           pasuv[3]=pasSav[3];
2739 #else 
2740           Standard_Real u1,v1,u2,v2;
2741           Standard_Real U1,V1,U2,V2;
2742           Standard_Integer nn=line->NbPoints();
2743           if(nn>2) { 
2744             line->Value(nn).Parameters(u1,v1,u2,v2);
2745             line->Value(nn-1).Parameters(U1,V1,U2,V2);
2746             pasuv[0]=Abs(u1-U1);
2747             pasuv[1]=Abs(v1-V1);
2748             pasuv[2]=Abs(u2-U2);
2749             pasuv[3]=Abs(v2-V2);
2750           }
2751 #endif
2752         }
2753         else Arrive = Standard_True;
2754     }
2755     else {
2756       pasuv[0]*=0.5;
2757       pasuv[1]*=0.5;
2758       pasuv[2]*=0.5;
2759       pasuv[3]*=0.5; 
2760     }
2761   }
2762 }
2763
2764 namespace {
2765   //OCC431(apo): modified ->
2766   static const Standard_Real CosRef2D =  Cos(M_PI/9.0),  AngRef2D = M_PI/2.0; 
2767
2768   static const Standard_Real d = 7.0;
2769 }
2770
2771 IntWalk_StatusDeflection  IntWalk_PWalking::TestDeflection()
2772
2773 // test if vector is observed by calculating an increase of vector 
2774 //     or the previous point and its tangent, the new calculated point and its  
2775 //     tangent; it is possible to find a cube passing by the 2 points and having as a 
2776 //     derivative the tangents of the intersection
2777 //     calculate the point with parameter 0.5 on cube=p1 
2778 //     calculate the medium point of 2 points of intersection=p2
2779 //   if arrow/2<=||p1p2||<= arrow consider that the vector is observed
2780 //   otherwise adjust the step depending on the ratio ||p1p2||/vector
2781 //   and the previous step 
2782 // test if in  2 tangent planes of surfaces there is no too great angle2d 
2783 // grand : if yes divide the step
2784 // test if there is no change of side
2785 //  
2786 {
2787   if(line->NbPoints() ==1 ) { 
2788     STATIC_BLOCAGE_SUR_PAS_TROP_GRAND=STATIC_PRECEDENT_INFLEXION=0;
2789   }
2790
2791   IntWalk_StatusDeflection Status = IntWalk_OK;
2792   Standard_Real FlecheCourante ,Ratio;
2793
2794
2795   const IntSurf_PntOn2S& CurrentPoint = myIntersectionOn2S.Point(); 
2796   //==================================================================================
2797   //=========               S t o p   o n   p o i n t                 ============
2798   //================================================================================== 
2799   if (myIntersectionOn2S.IsTangent())  { 
2800     return IntWalk_ArretSurPoint;  
2801   }
2802
2803   const gp_Dir& TgCourante = myIntersectionOn2S.Direction();
2804
2805   //==================================================================================
2806   //=========   R i s k   o f    i n f l e x i o n   p o i n t  ============
2807   //==================================================================================  
2808   if (TgCourante.Dot(previousd)<0) {
2809     //------------------------------------------------------------
2810     //-- Risk of inflexion point : Divide the step by 2
2811     //-- Initialize STATIC_PRECEDENT_INFLEXION so that 
2812     //-- at the next call to return Pas_OK if there is no 
2813     //-- more risk of the point of inflexion
2814     //------------------------------------------------------------
2815
2816     pasuv[0]*=0.5;
2817     pasuv[1]*=0.5;
2818     pasuv[2]*=0.5;
2819     pasuv[3]*=0.5;
2820     STATIC_PRECEDENT_INFLEXION+=3; 
2821     if (pasuv[0] < ResoU1 && pasuv[1] <ResoV1 && pasuv[2] <ResoU2 && pasuv[3] < ResoV2)
2822       return IntWalk_ArretSurPointPrecedent;
2823     else 
2824       return IntWalk_PasTropGrand;
2825   }
2826
2827   else {
2828     if(STATIC_PRECEDENT_INFLEXION  > 0) { 
2829       STATIC_PRECEDENT_INFLEXION -- ;
2830       return IntWalk_OK;
2831     }
2832   }
2833
2834   //==================================================================================
2835   //=========  D e t e c t    c o n f u s e d    P o in t s       ===========
2836   //==================================================================================
2837
2838   Standard_Real Dist = previousPoint.Value().
2839     SquareDistance(CurrentPoint.Value());
2840
2841
2842   if (Dist < tolconf*tolconf ) { 
2843     pasuv[0] = Max(5.*ResoU1,Min(1.5*pasuv[0],pasInit[0]));
2844     pasuv[1] = Max(5.*ResoV1,Min(1.5*pasuv[1],pasInit[1]));
2845     pasuv[2] = Max(5.*ResoU2,Min(1.5*pasuv[2],pasInit[2]));
2846     pasuv[3] = Max(5.*ResoV2,Min(1.5*pasuv[3],pasInit[3]));
2847     Status = IntWalk_PointConfondu;
2848   }
2849
2850   //==================================================================================
2851   Standard_Real Up1,Vp1,Uc1,Vc1,Du1,Dv1,AbsDu1,AbsDu2,AbsDv1,AbsDv2;
2852   Standard_Real Up2,Vp2,Uc2,Vc2,Du2,Dv2;
2853
2854   previousPoint.Parameters(Up1,Vp1,Up2,Vp2);
2855   CurrentPoint.Parameters(Uc1,Vc1,Uc2,Vc2);               
2856
2857   Du1 = Uc1 - Up1;   Dv1 = Vc1 - Vp1;
2858   Du2 = Uc2 - Up2;   Dv2 = Vc2 - Vp2;
2859
2860   AbsDu1 = Abs(Du1);
2861   AbsDu2 = Abs(Du2);
2862   AbsDv1 = Abs(Dv1);
2863   AbsDv2 = Abs(Dv2);
2864   //=================================================================================
2865   //====   S t e p   o f   p  r o g r e s s i o n (between previous and Current)   =======
2866   //=================================================================================
2867   if (   AbsDu1 < ResoU1 && AbsDv1 < ResoV1 
2868     && AbsDu2 < ResoU2 && AbsDv2 < ResoV2) {
2869       pasuv[0] = ResoU1; pasuv[1] = ResoV1; pasuv[2] = ResoU2; pasuv[3] = ResoV2;
2870       return(IntWalk_ArretSurPointPrecedent);
2871   }
2872   //==================================================================================
2873
2874   Standard_Real tolArea = 100.0;
2875   if (ResoU1 < Precision::PConfusion() ||
2876     ResoV1 < Precision::PConfusion() ||
2877     ResoU2 < Precision::PConfusion() ||
2878     ResoV2 < Precision::PConfusion() )
2879     tolArea =  tolArea*2.0;
2880
2881   Standard_Real Cosi1, CosRef1, Ang1, AngRef1, ResoUV1, Duv1, d1, tolCoeff1;   
2882   Standard_Real Cosi2, CosRef2, Ang2, AngRef2, ResoUV2, Duv2, d2, tolCoeff2;   
2883   Cosi1 = Du1*previousd1.X() + Dv1*previousd1.Y();
2884   Cosi2 = Du2*previousd2.X() + Dv2*previousd2.Y();
2885   Duv1 = Du1*Du1 + Dv1*Dv1;
2886   Duv2 = Du2*Du2 + Dv2*Dv2;
2887   ResoUV1 = ResoU1*ResoU1 + ResoV1*ResoV1;
2888   ResoUV2 = ResoU2*ResoU2 + ResoV2*ResoV2;
2889   //
2890   //modified by NIZNHY-PKV Wed Nov 13 12:25:44 2002 f
2891   //
2892   Standard_Real aMinDiv2=Precision::Confusion();
2893   aMinDiv2=aMinDiv2*aMinDiv2;
2894   //
2895   d1=d;
2896   if (Duv1>aMinDiv2)  {
2897     d1 = Abs(ResoUV1/Duv1);
2898     d1 = Min(Sqrt(d1)*tolArea, d);  
2899   } 
2900   //d1 = Abs(ResoUV1/Duv1); 
2901   //d1 = Min(Sqrt(d1)*tolArea,d);  
2902   //modified by NIZNHY-PKV Wed Nov 13 12:34:30 2002 t
2903   tolCoeff1 = Exp(d1);
2904   //
2905   //modified by NIZNHY-PKV Wed Nov 13 12:34:43 2002 f
2906   d2=d;
2907   if (Duv2>aMinDiv2) {
2908     d2 = Abs(ResoUV2/Duv2); 
2909     d2 = Min(Sqrt(d2)*tolArea,d); 
2910   }
2911   //d2 = Abs(ResoUV2/Duv2); 
2912   //d2 = Min(Sqrt(d2)*tolArea,d);  
2913   //modified by NIZNHY-PKV Wed Nov 13 12:34:53 2002 t
2914   tolCoeff2 = Exp(d2);
2915   CosRef1 = CosRef2D/tolCoeff1;
2916   CosRef2 = CosRef2D/tolCoeff2;
2917   //
2918   //==================================================================================
2919   //== The points are not confused :                                           ==
2920   //== 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 ==
2921   //==                           N o t    T o o    G r e a t (angle in space UV)    ==
2922   //==                           C h a n g e    o f    s i d e                ==
2923   //==================================================================================
2924   if (Status != IntWalk_PointConfondu) { 
2925     if(Cosi1*Cosi1 < CosRef1*Duv1 || Cosi2*Cosi2 < CosRef2*Duv2) {
2926       pasuv[0]*=0.5;  pasuv[1]*=0.5;  pasuv[2]*=0.5;  pasuv[3]*=0.5;
2927       if (pasuv[0]<ResoU1 && pasuv[1]<ResoV1 && pasuv[2]<ResoU2 && pasuv[3]<ResoV2) { 
2928         return(IntWalk_ArretSurPointPrecedent);
2929       }
2930       else {
2931         pasuv[0]*=0.5; pasuv[1]*=0.5; pasuv[2]*=0.5; pasuv[3]*=0.5;
2932         return(IntWalk_PasTropGrand);
2933       }
2934     }
2935     const gp_Dir2d& Tg2dcourante1 = myIntersectionOn2S.DirectionOnS1();
2936     const gp_Dir2d& Tg2dcourante2 = myIntersectionOn2S.DirectionOnS2();
2937     Cosi1 = Du1*Tg2dcourante1.X() + Dv1*Tg2dcourante1.Y();
2938     Cosi2 = Du2*Tg2dcourante2.X() + Dv2*Tg2dcourante2.Y();
2939     Ang1 = Abs(previousd1.Angle(Tg2dcourante1));  
2940     Ang2 = Abs(previousd2.Angle(Tg2dcourante2));  
2941     AngRef1 = AngRef2D*tolCoeff1;
2942     AngRef2 = AngRef2D*tolCoeff2;
2943     //-------------------------------------------------------
2944     //-- Test : Angle too great in space UV       -----
2945     //--        Change of  side                      -----
2946     //-------------------------------------------------------
2947     if(Cosi1*Cosi1 < CosRef1*Duv1 || Cosi2*Cosi2 < CosRef2*Duv2 || Ang1 > AngRef1 || Ang2 > AngRef2) {
2948       pasuv[0]*=0.5;  pasuv[1]*=0.5;  pasuv[2]*=0.5;  pasuv[3]*=0.5;
2949       if (pasuv[0]<ResoU1 && pasuv[1]<ResoV1 && pasuv[2]<ResoU2 && pasuv[3]<ResoV2) 
2950         return(IntWalk_ArretSurPoint);
2951       else 
2952         return(IntWalk_PasTropGrand);
2953     }
2954   }
2955   //<-OCC431(apo)
2956   //==================================================================================
2957   //== D e t e c t i o n   o f    :  Step Too Small 
2958   //==                               STEP TOO Great 
2959   //==================================================================================
2960
2961   //---------------------------------------
2962   //-- Estimate of the vector           --
2963   //---------------------------------------
2964   FlecheCourante =
2965     Sqrt(Abs((previousd.XYZ()-TgCourante.XYZ()).SquareModulus()*Dist))/8.;
2966
2967   if ( FlecheCourante<= fleche*0.5) {     //-- Current step too small
2968     if(FlecheCourante>1e-16) { 
2969       Ratio = 0.5*(fleche/FlecheCourante);
2970     }
2971     else { 
2972       Ratio = 10.0;
2973     }
2974     Standard_Real pasSu1 = pasuv[0];
2975     Standard_Real pasSv1 = pasuv[1];
2976     Standard_Real pasSu2 = pasuv[2];
2977     Standard_Real pasSv2 = pasuv[3];
2978
2979     //-- In  case if 
2980     //-- a point at U+DeltaU is required, ....
2981     //-- return a point at U + Epsilon
2982     //-- Epsilon << DeltaU.
2983
2984     if(pasuv[0]< AbsDu1) pasuv[0] = AbsDu1;
2985     if(pasuv[1]< AbsDv1) pasuv[1] = AbsDv1;
2986     if(pasuv[2]< AbsDu2) pasuv[2] = AbsDu2;
2987     if(pasuv[3]< AbsDv2) pasuv[3] = AbsDv2;
2988
2989     if(pasuv[0]<ResoU1) pasuv[0]=ResoU1;
2990     if(pasuv[1]<ResoV1) pasuv[1]=ResoV1;
2991     if(pasuv[2]<ResoU2) pasuv[2]=ResoU2;
2992     if(pasuv[3]<ResoV2) pasuv[3]=ResoV2;
2993     //-- if(Ratio>10.0 ) { Ratio=10.0; } 
2994     Standard_Real R1,R = pasInit[0]/pasuv[0];
2995     R1= pasInit[1]/pasuv[1];     if(R1<R) R=R1;
2996     R1= pasInit[2]/pasuv[2];     if(R1<R) R=R1;
2997     R1= pasInit[3]/pasuv[3];     if(R1<R) R=R1;
2998     if(Ratio > R) Ratio=R;
2999     pasuv[0] = Min(Ratio*pasuv[0],pasInit[0]);
3000     pasuv[1] = Min(Ratio*pasuv[1],pasInit[1]);
3001     pasuv[2] = Min(Ratio*pasuv[2],pasInit[2]);
3002     pasuv[3] = Min(Ratio*pasuv[3],pasInit[3]);
3003     if (pasuv[0] != pasSu1 || pasuv[2] != pasSu2|| 
3004       pasuv[1] != pasSv1 || pasuv[3] != pasSv2) {
3005         if(++STATIC_BLOCAGE_SUR_PAS_TROP_GRAND > 5) {
3006           STATIC_BLOCAGE_SUR_PAS_TROP_GRAND = 0;
3007           return IntWalk_PasTropGrand; 
3008         }
3009     }
3010     if(Status == IntWalk_OK) { 
3011       STATIC_BLOCAGE_SUR_PAS_TROP_GRAND=0;
3012       //-- Try to increase the step
3013     }
3014     return Status;
3015   }
3016   else {                                //-- CurrentVector > vector*0.5 
3017     if (FlecheCourante > fleche) {      //-- Current step too Great
3018       Ratio = fleche/FlecheCourante; 
3019       pasuv[0] = Ratio*pasuv[0];
3020       pasuv[1] = Ratio*pasuv[1];
3021       pasuv[2] = Ratio*pasuv[2];
3022       pasuv[3] = Ratio*pasuv[3];
3023       //if(++STATIC_BLOCAGE_SUR_PAS_TROP_GRAND > 5) {
3024       //        STATIC_BLOCAGE_SUR_PAS_TROP_GRAND = 0;
3025       return IntWalk_PasTropGrand; 
3026       //}
3027     }
3028     else {                             //-- vector/2  <  CurrentVector <= vector   
3029       Ratio = 0.75 * (fleche / FlecheCourante);
3030     }
3031   }
3032   pasuv[0] = Max(5.*ResoU1,Min(Min(Ratio*AbsDu1,pasuv[0]),pasInit[0]));
3033   pasuv[1] = Max(5.*ResoV1,Min(Min(Ratio*AbsDv1,pasuv[1]),pasInit[1]));
3034   pasuv[2] = Max(5.*ResoU2,Min(Min(Ratio*AbsDu2,pasuv[2]),pasInit[2]));
3035   pasuv[3] = Max(5.*ResoV2,Min(Min(Ratio*AbsDv2,pasuv[3]),pasInit[3]));
3036   if(Status == IntWalk_OK) STATIC_BLOCAGE_SUR_PAS_TROP_GRAND=0;
3037   return Status;
3038 }
3039
3040 Standard_Boolean IntWalk_PWalking::
3041 TestArret(const Standard_Boolean DejaReparti,
3042           TColStd_Array1OfReal& Param,
3043           IntImp_ConstIsoparametric&  ChoixIso)
3044
3045           //
3046           // test if the point of intersection set by these parameters remains in the 
3047           // natural domain of each square.
3048           // if the point outpasses reframe to find the best iso (border)
3049           // that intersects easiest the other square
3050           // otherwise test if closed line is present  
3051           // 
3052 {
3053   Standard_Real Uvd[4],Uvf[4],Epsuv[4],Duv[4],Uvp[4],dv,dv2,ParC[4];
3054   Standard_Real DPc,DPb;
3055   Standard_Integer i = 0, k = 0;
3056   Epsuv[0] = ResoU1;
3057   Epsuv[1] = ResoV1;
3058   Epsuv[2] = ResoU2;
3059   Epsuv[3] = ResoV2;
3060   previousPoint.Parameters(Uvp[0],Uvp[1],Uvp[2],Uvp[3]);
3061
3062   Standard_Real SolParam[4];
3063   myIntersectionOn2S.Point().Parameters(SolParam[0],SolParam[1],SolParam[2],SolParam[3]);
3064
3065   Standard_Boolean Trouve = Standard_False;
3066
3067   Uvd[0]=Um1;   Uvf[0]=UM1;   Uvd[1]=Vm1;   Uvf[1]=VM1;
3068   Uvd[2]=Um2;   Uvf[2]=UM2;   Uvd[3]=Vm2;   Uvf[3]=VM2;
3069
3070   Standard_Integer im1;
3071   for ( i = 1,im1 = 0;i<=4;i++,im1++) {
3072     switch(i) { 
3073     case 1: k=2; break;
3074     case 2: k=1; break;
3075     case 3: k=4; break;
3076     case 4: k=3; break;
3077     }
3078     if (Param(i) < (Uvd[im1]-Epsuv[im1]) ||
3079       SolParam[im1] < (Uvd[im1]-Epsuv[im1]))     //--     Current -----  Bound Inf -----  Previous
3080     {
3081       Trouve    = Standard_True;                   //-- 
3082       DPc       = Uvp[im1]-Param(i);               //--     Previous  - Current
3083       DPb       = Uvp[im1]-Uvd[im1];               //--     Previous  - Bound Inf
3084       ParC[im1] = Uvd[im1];                        //--     ParamCorrige
3085       dv        = Param(k)-Uvp[k-1];               //--     Current   - Previous (other Direction)
3086       dv2       = dv*dv;         
3087       if(dv2>RealEpsilon()) {                       //--    Progress at the other Direction ?
3088         Duv[im1]  = DPc*DPb + dv2;
3089         Duv[im1]  = Duv[im1]*Duv[im1]/(DPc*DPc+dv2)/(DPb*DPb+dv2);
3090       }
3091       else {
3092         Duv[im1]=-1.0;                              //--    If no progress, do not change  
3093       }                                             //--    the choice of iso 
3094     }   
3095     else if (Param(i) > (Uvf[im1] + Epsuv[im1]) ||
3096       SolParam[im1] > (Uvf[im1] + Epsuv[im1]))//--    Previous -----  Bound Sup -----  Current
3097     {
3098       Trouve    = Standard_True;                    //-- 
3099       DPc       = Param(i)-Uvp[im1];                //--     Current   - Previous
3100       DPb       = Uvf[im1]-Uvp[im1];                //--     Bound Sup - Previous 
3101       ParC[im1] = Uvf[im1];                         //--     Param Corrige
3102       dv        = Param(k)-Uvp[k-1];                //--     Current   - Previous (other Direction)
3103       dv2       = dv*dv;
3104       if(dv2>RealEpsilon()) {                       //--     Progress in other Direction ?
3105         Duv[im1]  =  DPc*DPb + dv2;
3106         Duv[im1]  = Duv[im1]*Duv[im1]/(DPc*DPc+dv2)/(DPb*DPb+dv2);
3107       }
3108       else {
3109         Duv[im1]=-1.0;                              //--    If no progress, do not change 
3110       }                                             //--    the choice of iso 
3111     }
3112     else { 
3113       Duv[im1]= -1.;
3114       ParC[im1]=Param(i);
3115     }
3116   }
3117
3118   if (Trouve) {
3119     //--------------------------------------------------
3120     //-- One of Parameters u1,v1,u2,v2 is outside of  --
3121     //-- the natural limits.                          -- 
3122     //-- Find the best direction of                   -- 
3123     //-- progress and reframe the parameters.        --
3124     //--------------------------------------------------
3125     Standard_Real ddv = -1.0;
3126     k=-1;
3127     for (i=0;i<=3;i++) {
3128       Param(i+1) = ParC[i];
3129       if(Duv[i]>ddv) { 
3130         ddv = Duv[i];
3131         k=i;
3132       }
3133     }
3134     if(k!=-1) { 
3135       ChoixIso   = ChoixRef[k];
3136     }
3137     else { 
3138       if((ParC[0]<=Uvd[0]+Epsuv[0]) || (ParC[0]>=Uvf[0]-Epsuv[0])) {
3139         ChoixIso = IntImp_UIsoparametricOnCaro1;
3140       }
3141       else if((ParC[1]<=Uvd[1]+Epsuv[1]) || (ParC[1]>=Uvf[1]-Epsuv[1])) {
3142         ChoixIso = IntImp_VIsoparametricOnCaro1;
3143       }
3144       else if((ParC[2]<=Uvd[2]+Epsuv[2]) || (ParC[2]>=Uvf[2]-Epsuv[2])) {
3145         ChoixIso = IntImp_UIsoparametricOnCaro2;
3146       }
3147       else if((ParC[3]<=Uvd[3]+Epsuv[3]) || (ParC[3]>=Uvf[3]-Epsuv[3])) {
3148         ChoixIso = IntImp_VIsoparametricOnCaro2;
3149       }
3150     }
3151     close = Standard_False;
3152     return Standard_True;
3153   }
3154   else 
3155   {  
3156     if (!DejaReparti) { // find if line closed
3157
3158       Standard_Real u,v;
3159       const IntSurf_PntOn2S& POn2S1=line->Value(1);
3160       //On S1
3161       POn2S1.ParametersOnS1(u,v);
3162       gp_Pnt2d P1uvS1(u,v);
3163       previousPoint.ParametersOnS1(u,v);
3164       gp_Pnt2d PrevuvS1(u,v);
3165       myIntersectionOn2S.Point().ParametersOnS1(u,v);
3166       gp_Pnt2d myIntersuvS1(u,v);
3167       Standard_Boolean close2dS1 = (P1uvS1.XY()-PrevuvS1.XY())*
3168         (P1uvS1.XY()-myIntersuvS1.XY()) < 0.0;
3169       //On S2
3170       POn2S1.ParametersOnS2(u,v);
3171       gp_Pnt2d P1uvS2(u,v);
3172       previousPoint.ParametersOnS2(u,v);
3173       gp_Pnt2d PrevuvS2(u,v);
3174       myIntersectionOn2S.Point().ParametersOnS2(u,v);
3175       gp_Pnt2d myIntersuvS2(u,v);
3176       Standard_Boolean close2dS2 = (P1uvS2.XY()-PrevuvS2.XY())*
3177         (P1uvS2.XY()-myIntersuvS2.XY()) < 0.0;
3178
3179       close = close2dS1 && close2dS2;
3180       return close;
3181     }
3182     else return Standard_False;
3183   }
3184 }
3185