0024472: Wrong section curves
[occt.git] / src / IntWalk / IntWalk_PWalking_1.gxx
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
7 // under the terms of the GNU Lesser General Public 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 //-----------------------------
16 //--  IntWalk_PWalking_1.gxx
17 //-- 
18
19 #include <Precision.hxx>
20 #include <math_FunctionSetRoot.hxx>
21
22 #include <Geom_Surface.hxx>
23
24 //#define KELARG 20.0
25
26 //==================================================================================
27 // function : IntWalk_PWalking::IntWalk_PWalking
28 // purpose  :
29 // estimate of max step : To avoid abrupt changes 
30 // during change of isos 
31 //==================================================================================
32 void ComputePasInit(Standard_Real *pasuv,
33                     Standard_Real Um1,Standard_Real UM1,
34                     Standard_Real Vm1,Standard_Real VM1,
35                     Standard_Real Um2,Standard_Real UM2,
36                     Standard_Real Vm2,Standard_Real VM2,
37                     Standard_Real _Um1,Standard_Real _UM1,
38                     Standard_Real _Vm1,Standard_Real _VM1,
39                     Standard_Real _Um2,Standard_Real _UM2,
40                     Standard_Real _Vm2,Standard_Real _VM2,
41                     const ThePSurface& ,
42                     const ThePSurface& ,
43                     const Standard_Real Increment) 
44
45   Standard_Real du1=Abs(UM1-Um1);
46   Standard_Real dv1=Abs(VM1-Vm1);
47   Standard_Real du2=Abs(UM2-Um2);
48   Standard_Real dv2=Abs(VM2-Vm2);
49
50   Standard_Real _du1=Abs(_UM1-_Um1);
51   Standard_Real _dv1=Abs(_VM1-_Vm1);
52   Standard_Real _du2=Abs(_UM2-_Um2);
53   Standard_Real _dv2=Abs(_VM2-_Vm2);
54
55   //-- limit the reduction of uv box estimate to 0.01 natural box
56   //--  du1 : On box of Inter
57   //-- _du1 : On parametric space
58   if(_du1<1e50 && du1<0.01*_du1) du1=0.01*_du1;
59   if(_dv1<1e50 && dv1<0.01*_dv1) dv1=0.01*_dv1;
60   if(_du2<1e50 && du2<0.01*_du2) du2=0.01*_du2;
61   if(_dv2<1e50 && dv2<0.01*_dv2) dv2=0.01*_dv2;
62
63   pasuv[0]=Increment*du1;
64   pasuv[1]=Increment*dv1;
65   pasuv[2]=Increment*du2;
66   pasuv[3]=Increment*dv2;
67 }
68 //==================================================================================
69 // function : IntWalk_PWalking::IntWalk_PWalking
70 // purpose  : 
71 //==================================================================================
72 IntWalk_PWalking::IntWalk_PWalking(const ThePSurface& Caro1,
73                                    const ThePSurface& Caro2,
74                                    const Standard_Real TolTangency,
75                                    const Standard_Real Epsilon,
76                                    const Standard_Real Deflection,
77                                    const Standard_Real Increment ) 
78                                    :
79
80 done(Standard_True),
81 close(Standard_False),
82 fleche(Deflection),
83 tolconf(Epsilon),
84 sensCheminement(1),
85 myIntersectionOn2S(Caro1,Caro2,TolTangency),
86 STATIC_BLOCAGE_SUR_PAS_TROP_GRAND(0),
87 STATIC_PRECEDENT_INFLEXION(0)
88 {
89   Standard_Real KELARG=20.;
90   //
91   pasMax=Increment*0.2; //-- June 25 99 after problems with precision 
92   Um1 = ThePSurfaceTool::FirstUParameter(Caro1);
93   Vm1 = ThePSurfaceTool::FirstVParameter(Caro1);
94   UM1 = ThePSurfaceTool::LastUParameter(Caro1);
95   VM1 = ThePSurfaceTool::LastVParameter(Caro1);
96
97   Um2 = ThePSurfaceTool::FirstUParameter(Caro2);
98   Vm2 = ThePSurfaceTool::FirstVParameter(Caro2);
99   UM2 = ThePSurfaceTool::LastUParameter(Caro2);
100   VM2 = ThePSurfaceTool::LastVParameter(Caro2);
101
102   ResoU1 = ThePSurfaceTool::UResolution(Caro1,Precision::Confusion());
103   ResoV1 = ThePSurfaceTool::VResolution(Caro1,Precision::Confusion());
104
105   ResoU2 = ThePSurfaceTool::UResolution(Caro2,Precision::Confusion());
106   ResoV2 = ThePSurfaceTool::VResolution(Caro2,Precision::Confusion());
107
108   Standard_Real NEWRESO;
109   Standard_Real MAXVAL;
110   Standard_Real MAXVAL2;
111   //
112   MAXVAL  = Abs(Um1);  MAXVAL2 = Abs(UM1);
113   if(MAXVAL2 > MAXVAL) MAXVAL = MAXVAL2;
114   NEWRESO = ResoU1 * MAXVAL ;
115   if(NEWRESO > ResoU1 &&NEWRESO<10) {    ResoU1 = NEWRESO;  }
116
117
118   MAXVAL  = Abs(Um2);   MAXVAL2 = Abs(UM2);
119   if(MAXVAL2 > MAXVAL) MAXVAL = MAXVAL2;
120   NEWRESO = ResoU2 * MAXVAL ;
121   if(NEWRESO > ResoU2 && NEWRESO<10) {     ResoU2 = NEWRESO;  }
122
123
124   MAXVAL  = Abs(Vm1);  MAXVAL2 = Abs(VM1);
125   if(MAXVAL2 > MAXVAL) MAXVAL = MAXVAL2;
126   NEWRESO = ResoV1 * MAXVAL ;
127   if(NEWRESO > ResoV1 && NEWRESO<10) {     ResoV1 = NEWRESO;  }
128
129
130   MAXVAL  = Abs(Vm2);  MAXVAL2 = Abs(VM2);
131   if(MAXVAL2 > MAXVAL) MAXVAL = MAXVAL2;
132   NEWRESO = ResoV2 * MAXVAL ;
133   if(NEWRESO > ResoV2 && NEWRESO<10) {     ResoV2 = NEWRESO;  }
134
135   pasuv[0]=pasMax*Abs(UM1-Um1);
136   pasuv[1]=pasMax*Abs(VM1-Vm1);
137   pasuv[2]=pasMax*Abs(UM2-Um2);
138   pasuv[3]=pasMax*Abs(VM2-Vm2);
139
140   if(ResoU1>0.0001*pasuv[0]) ResoU1=0.00001*pasuv[0];
141   if(ResoV1>0.0001*pasuv[1]) ResoV1=0.00001*pasuv[1];
142   if(ResoU2>0.0001*pasuv[2]) ResoU2=0.00001*pasuv[2];
143   if(ResoV2>0.0001*pasuv[3]) ResoV2=0.00001*pasuv[3];
144
145
146   if(ThePSurfaceTool::IsUPeriodic(Caro1)==Standard_False) { 
147     //UM1+=KELARG*pasuv[0];  Um1-=KELARG*pasuv[0];
148   }
149   else { 
150     Standard_Real t = UM1-Um1; 
151     if(t<ThePSurfaceTool::UPeriod(Caro1)) { 
152       t=0.5*(ThePSurfaceTool::UPeriod(Caro1)-t);
153       t=(t>KELARG*pasuv[0])? KELARG*pasuv[0] : t;
154       UM1+=t;  Um1-=t;
155     }
156   }
157
158   if(ThePSurfaceTool::IsVPeriodic(Caro1)==Standard_False) { 
159     //VM1+=KELARG*pasuv[1];  Vm1-=KELARG*pasuv[1];
160   }
161   else { 
162     Standard_Real t = VM1-Vm1; 
163     if(t<ThePSurfaceTool::VPeriod(Caro1)) { 
164       t=0.5*(ThePSurfaceTool::VPeriod(Caro1)-t);
165       t=(t>KELARG*pasuv[1])? KELARG*pasuv[1] : t;
166       VM1+=t;  Vm1-=t;
167     }
168   }
169
170   if(ThePSurfaceTool::IsUPeriodic(Caro2)==Standard_False) { 
171     //UM2+=KELARG*pasuv[2];  Um2-=KELARG*pasuv[2];
172   }
173   else { 
174     Standard_Real t = UM2-Um2; 
175     if(t<ThePSurfaceTool::UPeriod(Caro2)) { 
176       t=0.5*(ThePSurfaceTool::UPeriod(Caro2)-t);
177       t=(t>KELARG*pasuv[2])? KELARG*pasuv[2] : t;
178       UM2+=t;  Um2-=t;
179     }
180   }
181
182   if(ThePSurfaceTool::IsVPeriodic(Caro2)==Standard_False) {   
183     //VM2+=KELARG*pasuv[3];  Vm2-=KELARG*pasuv[3];
184   }
185   else { 
186     Standard_Real t = VM2-Vm2; 
187     if(t<ThePSurfaceTool::VPeriod(Caro2)) { 
188       t=0.5*(ThePSurfaceTool::VPeriod(Caro2)-t);
189       t=(t>KELARG*pasuv[3])? KELARG*pasuv[3] : t;
190       VM2+=t;  Vm2-=t;
191     }
192   }
193
194   //-- ComputePasInit(pasuv,Um1,UM1,Vm1,VM1,Um2,UM2,Vm2,VM2,Caro1,Caro2);
195
196   for (Standard_Integer i = 0; i<=3;i++) {
197     if(pasuv[i]>10) 
198       pasuv[i] = 10; 
199     pasInit[i] = pasSav[i] = pasuv[i]; 
200   }
201
202
203 }
204 //==================================================================================
205 // function : IntWalk_PWalking
206 // purpose  : 
207 //==================================================================================
208 IntWalk_PWalking::IntWalk_PWalking(const ThePSurface& Caro1,
209                                    const ThePSurface& Caro2,
210                                    const Standard_Real TolTangency,
211                                    const Standard_Real Epsilon,
212                                    const Standard_Real Deflection,
213                                    const Standard_Real Increment, 
214                                    const Standard_Real U1,
215                                    const Standard_Real V1,
216                                    const Standard_Real U2, 
217                                    const Standard_Real V2)
218                                    :
219
220 done(Standard_True),
221 close(Standard_False),
222 fleche(Deflection),
223 tolconf(Epsilon),
224 sensCheminement(1),       
225 myIntersectionOn2S(Caro1,Caro2,TolTangency),
226 STATIC_BLOCAGE_SUR_PAS_TROP_GRAND(0),
227 STATIC_PRECEDENT_INFLEXION(0)
228 {
229   Standard_Real KELARG=20.;
230   //
231   pasMax=Increment*0.2; //-- June 25 99 after problems with precision 
232   //
233   Um1 = ThePSurfaceTool::FirstUParameter(Caro1);
234   Vm1 = ThePSurfaceTool::FirstVParameter(Caro1);
235   UM1 = ThePSurfaceTool::LastUParameter(Caro1);
236   VM1 = ThePSurfaceTool::LastVParameter(Caro1);
237
238   Um2 = ThePSurfaceTool::FirstUParameter(Caro2);
239   Vm2 = ThePSurfaceTool::FirstVParameter(Caro2);
240   UM2 = ThePSurfaceTool::LastUParameter(Caro2);
241   VM2 = ThePSurfaceTool::LastVParameter(Caro2);
242
243   ResoU1 = ThePSurfaceTool::UResolution(Caro1,Precision::Confusion());
244   ResoV1 = ThePSurfaceTool::VResolution(Caro1,Precision::Confusion());
245
246   ResoU2 = ThePSurfaceTool::UResolution(Caro2,Precision::Confusion());
247   ResoV2 = ThePSurfaceTool::VResolution(Caro2,Precision::Confusion());
248   //
249   Standard_Real NEWRESO, MAXVAL, MAXVAL2;
250   //
251   MAXVAL  = Abs(Um1);  
252   MAXVAL2 = Abs(UM1);
253   if(MAXVAL2 > MAXVAL) {
254     MAXVAL = MAXVAL2;
255   }
256   NEWRESO = ResoU1 * MAXVAL ;
257   if(NEWRESO > ResoU1) {
258     ResoU1 = NEWRESO;  
259   }
260   //
261   MAXVAL  = Abs(Um2);   
262   MAXVAL2 = Abs(UM2);
263   if(MAXVAL2 > MAXVAL){
264     MAXVAL = MAXVAL2;
265   }  
266   NEWRESO = ResoU2 * MAXVAL ;
267   if(NEWRESO > ResoU2) {
268     ResoU2 = NEWRESO;  
269   }
270   //
271   MAXVAL  = Abs(Vm1);  
272   MAXVAL2 = Abs(VM1);
273   if(MAXVAL2 > MAXVAL) {
274     MAXVAL = MAXVAL2;
275   }
276   NEWRESO = ResoV1 * MAXVAL ;
277   if(NEWRESO > ResoV1) {    
278     ResoV1 = NEWRESO; 
279   }
280   //
281   MAXVAL  = Abs(Vm2);  
282   MAXVAL2 = Abs(VM2);
283   if(MAXVAL2 > MAXVAL){
284     MAXVAL = MAXVAL2;
285   }  
286   NEWRESO = ResoV2 * MAXVAL ;
287   if(NEWRESO > ResoV2) {  
288     ResoV2 = NEWRESO;
289   }
290   //
291   pasuv[0]=pasMax*Abs(UM1-Um1);
292   pasuv[1]=pasMax*Abs(VM1-Vm1);
293   pasuv[2]=pasMax*Abs(UM2-Um2);
294   pasuv[3]=pasMax*Abs(VM2-Vm2);
295   //
296   if(ThePSurfaceTool::IsUPeriodic(Caro1)==Standard_False) { 
297     UM1+=KELARG*pasuv[0];  
298     Um1-=KELARG*pasuv[0];
299   }
300   else { 
301     Standard_Real t = UM1-Um1; 
302     if(t<ThePSurfaceTool::UPeriod(Caro1)) { 
303       t=0.5*(ThePSurfaceTool::UPeriod(Caro1)-t);
304       t=(t>KELARG*pasuv[0])? KELARG*pasuv[0] : t;
305       UM1+=t;  
306       Um1-=t;
307     }
308   }
309   //
310   if(ThePSurfaceTool::IsVPeriodic(Caro1)==Standard_False) { 
311     VM1+=KELARG*pasuv[1];
312     Vm1-=KELARG*pasuv[1];
313   }
314   else { 
315     Standard_Real t = VM1-Vm1; 
316     if(t<ThePSurfaceTool::VPeriod(Caro1)) { 
317       t=0.5*(ThePSurfaceTool::VPeriod(Caro1)-t);
318       t=(t>KELARG*pasuv[1])? KELARG*pasuv[1] : t;
319       VM1+=t;  Vm1-=t;
320     }
321   }
322   //
323   if(ThePSurfaceTool::IsUPeriodic(Caro2)==Standard_False) { 
324     UM2+=KELARG*pasuv[2];  
325     Um2-=KELARG*pasuv[2];
326   }
327   else { 
328     Standard_Real t = UM2-Um2; 
329     if(t<ThePSurfaceTool::UPeriod(Caro2)) { 
330       t=0.5*(ThePSurfaceTool::UPeriod(Caro2)-t);
331       t=(t>KELARG*pasuv[2])? KELARG*pasuv[2] : t;
332       UM2+=t;  
333       Um2-=t;
334     }
335   }
336
337   if(ThePSurfaceTool::IsVPeriodic(Caro2)==Standard_False) {   
338     VM2+=KELARG*pasuv[3];  
339     Vm2-=KELARG*pasuv[3];
340   }
341   else { 
342     Standard_Real t = VM2-Vm2; 
343     if(t<ThePSurfaceTool::UPeriod(Caro2)) { 
344       t=0.5*(ThePSurfaceTool::VPeriod(Caro2)-t);
345       t=(t>KELARG*pasuv[3])? KELARG*pasuv[3] : t;
346       VM2+=t;  
347       Vm2-=t;
348     }
349   }
350   //-- ComputePasInit(pasuv,Um1,UM1,Vm1,VM1,Um2,UM2,Vm2,VM2,Caro1,Caro2);
351
352   for (Standard_Integer i = 0; i<=3;i++) {
353     pasInit[i] = pasSav[i] = pasuv[i]; 
354   }  
355
356   if(ResoU1>0.0001*pasuv[0]) ResoU1=0.00001*pasuv[0];
357   if(ResoV1>0.0001*pasuv[1]) ResoV1=0.00001*pasuv[1];
358   if(ResoU2>0.0001*pasuv[2]) ResoU2=0.00001*pasuv[2];
359   if(ResoV2>0.0001*pasuv[3]) ResoV2=0.00001*pasuv[3];
360   //
361   TColStd_Array1OfReal Par(1,4);
362   Par(1) = U1;
363   Par(2) = V1;
364   Par(3) = U2;
365   Par(4) = V2;
366   Perform(Par);
367 }
368
369 //==================================================================================
370 // function : PerformFirstPoint
371 // purpose  : 
372 //==================================================================================
373 Standard_Boolean IntWalk_PWalking::PerformFirstPoint  (const TColStd_Array1OfReal& ParDep,
374                                                        IntSurf_PntOn2S& FirstPoint)   
375 {
376   sensCheminement = 1;
377   close = Standard_False;
378   //
379   Standard_Integer i;
380   TColStd_Array1OfReal Param(1,4);
381   //
382   for (i=1; i<=4; ++i) {
383     Param(i) = ParDep(i);
384   }
385   //-- calculate the first solution point
386   math_FunctionSetRoot  Rsnld(myIntersectionOn2S.Function());
387   //
388   myIntersectionOn2S.Perform(Param,Rsnld);
389   if (!myIntersectionOn2S.IsDone())  { 
390     return Standard_False;
391   }
392
393   if (myIntersectionOn2S.IsEmpty()) {
394     return Standard_False;
395   }
396
397   FirstPoint = myIntersectionOn2S.Point();
398   return Standard_True;
399 }
400 //==================================================================================
401 // function : Perform
402 // purpose  : 
403 //==================================================================================
404 void IntWalk_PWalking::Perform(const TColStd_Array1OfReal& ParDep)    
405 {
406   Perform(ParDep,Um1,Vm1,Um2,Vm2,UM1,VM1,UM2,VM2);
407 }
408 //==================================================================================
409 // function : Perform
410 // purpose  : 
411 //==================================================================================
412 void IntWalk_PWalking::Perform(const TColStd_Array1OfReal& ParDep,
413                                const Standard_Real u1min,
414                                const Standard_Real v1min,
415                                const Standard_Real u2min,
416                                const Standard_Real v2min,
417                                const Standard_Real u1max,
418                                const Standard_Real v1max,
419                                const Standard_Real u2max,
420                                const Standard_Real v2max)
421 {
422   const Standard_Real aSQDistMax = 1.0e-14;
423   //xf
424
425   Standard_Integer NbPasOKConseq=0;
426   Standard_Real pasMaxSV[4], aTmp;
427   TColStd_Array1OfReal Param(1,4);
428   IntImp_ConstIsoparametric ChoixIso;
429   //xt
430   //
431   done = Standard_False;
432   //
433   // Caro1 and Caro2
434   const ThePSurface& Caro1 =myIntersectionOn2S.Function().AuxillarSurface1();
435   const ThePSurface& Caro2 =myIntersectionOn2S.Function().AuxillarSurface2();
436   //
437   const Standard_Real UFirst1 = ThePSurfaceTool::FirstUParameter(Caro1);
438   const Standard_Real VFirst1 = ThePSurfaceTool::FirstVParameter(Caro1);
439   const Standard_Real ULast1  = ThePSurfaceTool::LastUParameter (Caro1);
440   const Standard_Real VLast1  = ThePSurfaceTool::LastVParameter (Caro1);
441
442   const Standard_Real UFirst2 = ThePSurfaceTool::FirstUParameter(Caro2);
443   const Standard_Real VFirst2 = ThePSurfaceTool::FirstVParameter(Caro2);
444   const Standard_Real ULast2  = ThePSurfaceTool::LastUParameter (Caro2);
445   const Standard_Real VLast2  = ThePSurfaceTool::LastVParameter (Caro2);
446   //
447   ComputePasInit(pasuv,u1min,u1max,v1min,v1max,u2min,u2max,v2min,v2max,
448     Um1,UM1,Vm1,VM1,Um2,UM2,Vm2,VM2,Caro1,Caro2,pasMax+pasMax);
449   //
450   if(pasuv[0]<100.0*ResoU1) {
451     pasuv[0]=100.0*ResoU1; 
452   }
453   if(pasuv[1]<100.0*ResoV1) {
454     pasuv[1]=100.0*ResoV1; 
455   }
456   if(pasuv[2]<100.0*ResoU2) {
457     pasuv[2]=100.0*ResoU2;
458   }
459   if(pasuv[3]<100.0*ResoV2) {
460     pasuv[3]=100.0*ResoV2;
461   }
462   //
463   for (Standard_Integer i=0; i<4; ++i)
464   {
465     if(pasuv[i]>10)
466     {
467       pasuv[i] = 10;
468     }
469
470     pasInit[i] = pasSav[i] = pasuv[i]; 
471   }
472   //
473   line = new IntSurf_LineOn2S ();
474   //
475   for (Standard_Integer i=1; i<=4; ++i)
476   {
477     aTmp=ParDep(i);
478     Param(i)=ParDep(i);
479   }
480   //-- reproduce steps uv connected to surfaces Caro1 and Caro2
481   //-- pasuv[] and pasSav[] are modified during the marching
482   for(Standard_Integer i = 0; i < 4; ++i)
483   {
484     pasMaxSV[i] = pasSav[i] = pasuv[i] = pasInit[i]; 
485   }
486
487   //-- calculate the first solution point
488   math_FunctionSetRoot  Rsnld(myIntersectionOn2S.Function());
489   //
490   ChoixIso = myIntersectionOn2S.Perform(Param,Rsnld);
491   if (!myIntersectionOn2S.IsDone())
492   {
493     return;
494   }
495
496   //
497   if (myIntersectionOn2S.IsEmpty())
498   {
499     return;
500   }
501   //
502   if(myIntersectionOn2S.IsTangent())
503   {
504     return;
505   }
506   //
507   Standard_Boolean Arrive, DejaReparti;
508   const Standard_Integer RejectIndexMAX = 250000;
509   Standard_Integer IncKey, RejectIndex;
510   gp_Pnt pf,pl;
511   //
512   DejaReparti = Standard_False;
513   IncKey = 0;
514   RejectIndex = 0;
515   //
516   previousPoint = myIntersectionOn2S.Point();
517   previoustg = Standard_False;
518   previousd  = myIntersectionOn2S.Direction();
519   previousd1 = myIntersectionOn2S.DirectionOnS1();
520   previousd2 = myIntersectionOn2S.DirectionOnS2();
521   indextg = 1;
522   tgdir   = previousd;
523   firstd1 = previousd1;
524   firstd2 = previousd2;
525   tgfirst = tglast = Standard_False;
526   choixIsoSav  =  ChoixIso;
527   //------------------------------------------------------------
528   //-- Test if the first point of marching corresponds 
529   //-- to a point on borders. 
530   //-- In this case, DejaReparti is initialized as True
531   //-- 
532   pf = previousPoint.Value();
533   Standard_Boolean bTestFirstPoint = Standard_True;
534
535   previousPoint.Parameters(Param(1),Param(2),Param(3),Param(4));
536   AddAPoint(line,previousPoint);
537   //
538   IntWalk_StatusDeflection Status = IntWalk_OK;
539   Standard_Boolean NoTestDeflection = Standard_False;
540   Standard_Real SvParam[4], f;
541   Standard_Integer LevelOfEmptyInmyIntersectionOn2S=0;
542   Standard_Integer LevelOfPointConfondu = 0; 
543   Standard_Integer LevelOfIterWithoutAppend = -1;
544   //
545   Arrive = Standard_False;
546   while(!Arrive) //010
547   {
548     LevelOfIterWithoutAppend++;
549     if(LevelOfIterWithoutAppend>20)
550     {
551       Arrive = Standard_True; 
552       if(DejaReparti) {
553         break;
554       }
555       RepartirOuDiviser(DejaReparti,ChoixIso,Arrive);
556       LevelOfIterWithoutAppend = 0;
557     }
558     //
559     // compute f
560     f = 0.;
561     switch (ChoixIso) { 
562       case IntImp_UIsoparametricOnCaro1: f = Abs(previousd1.X()); break;
563       case IntImp_VIsoparametricOnCaro1: f = Abs(previousd1.Y()); break;
564       case IntImp_UIsoparametricOnCaro2: f = Abs(previousd2.X()); break;
565       case IntImp_VIsoparametricOnCaro2: f = Abs(previousd2.Y()); break;
566       default:break;
567     }
568     //
569     if(f<0.1) {
570       f=0.1;
571     }
572     //
573     previousPoint.Parameters(Param(1),Param(2),Param(3),Param(4));
574     //
575     //--ofv.begin
576     Standard_Real aIncKey, aEps, dP1, dP2, dP3, dP4;
577     //
578     dP1 = sensCheminement * pasuv[0] * previousd1.X() /f;
579     dP2 = sensCheminement * pasuv[1] * previousd1.Y() /f;
580     dP3 = sensCheminement * pasuv[2] * previousd2.X() /f; 
581     dP4 = sensCheminement * pasuv[3] * previousd2.Y() /f;
582     //
583     aIncKey=5.*(Standard_Real)IncKey;
584     aEps=1.e-7;
585     if(ChoixIso == IntImp_UIsoparametricOnCaro1 && Abs(dP1) < aEps)
586     {
587       dP1 *= aIncKey;
588     }
589
590     if(ChoixIso == IntImp_VIsoparametricOnCaro1 && Abs(dP2) < aEps)
591     {
592       dP2 *= aIncKey;
593     }
594
595     if(ChoixIso == IntImp_UIsoparametricOnCaro2 && Abs(dP3) < aEps)
596     {
597       dP3 *= aIncKey;
598     }
599
600     if(ChoixIso == IntImp_VIsoparametricOnCaro2 && Abs(dP4) < aEps)
601     {
602       dP4 *= aIncKey;
603     }
604     //--ofv.end
605     //
606     Param(1) += dP1;
607     Param(2) += dP2;
608     Param(3) += dP3; 
609     Param(4) += dP4;
610     //==========================
611     SvParam[0]=Param(1); 
612     SvParam[1]=Param(2);
613     SvParam[2]=Param(3);
614     SvParam[3]=Param(4);
615     //
616     ChoixIso= myIntersectionOn2S.Perform(Param, Rsnld, ChoixIso);                  
617     //
618     if (!myIntersectionOn2S.IsDone())
619     {
620       //end of line, division
621       Arrive = Standard_False;
622       Param(1)=SvParam[0]; 
623       Param(2)=SvParam[1]; 
624       Param(3)=SvParam[2];
625       Param(4)=SvParam[3];
626       RepartirOuDiviser(DejaReparti, ChoixIso, Arrive);
627     }
628     else  //009 
629     {
630       //== Calculation of exact point from Param(.) is possible
631       if (myIntersectionOn2S.IsEmpty())
632       {
633         Standard_Real u1,v1,u2,v2;
634         previousPoint.Parameters(u1,v1,u2,v2);
635         //
636         Arrive = Standard_False;
637         if(u1<UFirst1 || u1>ULast1)
638         {
639           Arrive=Standard_True;
640         }       
641
642         if(u2<UFirst2 || u2>ULast2)
643         {
644           Arrive=Standard_True;
645         }
646
647         if(v1<VFirst1 || v1>VLast1)
648         {
649           Arrive=Standard_True;
650         }
651
652         if(v2<VFirst2 || v2>VLast2)
653         {
654           Arrive=Standard_True;
655         }
656
657         RepartirOuDiviser(DejaReparti,ChoixIso,Arrive);
658         LevelOfEmptyInmyIntersectionOn2S++;
659         //
660         if(LevelOfEmptyInmyIntersectionOn2S>10)
661         {
662           pasuv[0]=pasSav[0]; 
663           pasuv[1]=pasSav[1]; 
664           pasuv[2]=pasSav[2]; 
665           pasuv[3]=pasSav[3];
666         }
667       }
668       else //008
669       {
670         //============================================================
671         //== A point has been found :  T E S T   D E F L E C T I O N 
672         //============================================================
673         if(NoTestDeflection)
674         {
675           NoTestDeflection = Standard_False;
676         }                 
677         else
678         {
679           if(--LevelOfEmptyInmyIntersectionOn2S<=0)
680           {
681             LevelOfEmptyInmyIntersectionOn2S=0;
682             if(LevelOfIterWithoutAppend < 10)
683             {
684               Status = TestDeflection();
685             }                   
686             else
687             {
688               pasuv[0]*=0.5; 
689               pasuv[1]*=0.5; 
690               pasuv[2]*=0.5; 
691               pasuv[3]*=0.5;
692             }
693           }
694         }
695
696         //============================================================
697         //==       T r a i t e m e n t   s u r   S t a t u s        ==
698         //============================================================
699         if(LevelOfPointConfondu > 5)
700         { 
701           Status = IntWalk_ArretSurPoint; 
702           LevelOfPointConfondu = 0;  
703         }
704         //
705         if(Status==IntWalk_OK)
706         { 
707           NbPasOKConseq++;
708           if(NbPasOKConseq >= 5)
709           {
710             NbPasOKConseq=0;
711             Standard_Boolean pastroppetit;
712             Standard_Real t;
713             //
714             do
715             {
716               pastroppetit=Standard_True;
717               //
718               if(pasuv[0]<pasInit[0])
719               {
720                 t = (pasInit[0]-pasuv[0])*0.25;
721                 if(t>0.1*pasInit[0])
722                 {
723                   t=0.1*pasuv[0];
724                 }
725
726                 pasuv[0]+=t; 
727                 pastroppetit=Standard_False;
728               }
729
730               if(pasuv[1]<pasInit[1])
731               {
732                 t = (pasInit[1]-pasuv[1])*0.25;
733                 if(t>0.1*pasInit[1]) {
734                   t=0.1*pasuv[1];
735                 }               
736
737                 pasuv[1]+=t; 
738                 pastroppetit=Standard_False;
739               }
740
741               if(pasuv[2]<pasInit[2])
742               {
743                 t = (pasInit[2]-pasuv[2])*0.25;
744                 if(t>0.1*pasInit[2])
745                 {
746                   t=0.1*pasuv[2];
747                 }
748
749                 pasuv[2]+=t; 
750                 pastroppetit=Standard_False;
751               }
752
753               if(pasuv[3]<pasInit[3])
754               {
755                 t = (pasInit[3]-pasuv[3])*0.25;
756                 if(t>0.1*pasInit[3]) {
757                   t=0.1*pasuv[3];
758                 }
759                 pasuv[3]+=t; 
760                 pastroppetit=Standard_False;
761               }
762               if(pastroppetit)
763               {
764                 if(pasMax<0.1)
765                 {
766                   pasMax*=1.1;
767                   pasInit[0]*=1.1; 
768                   pasInit[1]*=1.1; 
769                   pasInit[2]*=1.1; 
770                   pasInit[3]*=1.1; 
771                 }
772                 else
773                 {
774                   pastroppetit=Standard_False;
775                 }
776               }
777             }
778             while(pastroppetit);
779           }
780         }//Status==IntWalk_OK
781         else
782           NbPasOKConseq=0;
783
784         //
785         switch(Status)//007 
786         {
787         case IntWalk_ArretSurPointPrecedent:
788           {
789             Arrive = Standard_False;
790             RepartirOuDiviser(DejaReparti, ChoixIso, Arrive);
791             break;
792           }
793         case IntWalk_PasTropGrand:
794           {
795             Param(1)=SvParam[0];
796             Param(2)=SvParam[1]; 
797             Param(3)=SvParam[2]; 
798             Param(4)=SvParam[3];
799
800             if(LevelOfIterWithoutAppend > 5)
801             {
802               if(pasSav[0]<pasInit[0])
803               {
804                 pasInit[0]-=(pasInit[0]-pasSav[0])*0.25;
805                 LevelOfIterWithoutAppend=0;
806               }
807
808               if(pasSav[1]<pasInit[1])
809               {
810                 pasInit[1]-=(pasInit[1]-pasSav[1])*0.25;
811                 LevelOfIterWithoutAppend=0;
812               }
813
814               if(pasSav[2]<pasInit[2])
815               {
816                 pasInit[2]-=(pasInit[2]-pasSav[2])*0.25;
817                 LevelOfIterWithoutAppend=0;
818               }
819
820               if(pasSav[3]<pasInit[3])
821               {
822                 pasInit[3]-=(pasInit[3]-pasSav[3])*0.25;
823                 LevelOfIterWithoutAppend=0;
824               }
825             }
826
827             break;
828           }
829         case IntWalk_PointConfondu:
830           {
831             LevelOfPointConfondu++;
832
833             if(LevelOfPointConfondu>5)
834             {
835               Standard_Boolean pastroppetit;
836               //
837               do
838               {
839                 pastroppetit=Standard_True;
840
841                 if(pasuv[0]<pasInit[0])
842                 {
843                   pasuv[0]+=(pasInit[0]-pasuv[0])*0.25;
844                   pastroppetit=Standard_False;
845                 }
846
847                 if(pasuv[1]<pasInit[1])
848                 {
849                   pasuv[1]+=(pasInit[1]-pasuv[1])*0.25;
850                   pastroppetit=Standard_False;
851                 }
852
853                 if(pasuv[2]<pasInit[2])
854                 {
855                   pasuv[2]+=(pasInit[2]-pasuv[2])*0.25;
856                   pastroppetit=Standard_False; 
857                 }
858
859                 if(pasuv[3]<pasInit[3])
860                 {
861                   pasuv[3]+=(pasInit[3]-pasuv[3])*0.25;
862                   pastroppetit=Standard_False;
863                 }
864
865                 if(pastroppetit)
866                 {
867                   if(pasMax<0.1)
868                   {
869                     pasMax*=1.1;
870                     pasInit[0]*=1.1;
871                     pasInit[1]*=1.1;
872                     pasInit[2]*=1.1;
873                     pasInit[3]*=1.1; 
874                   }
875                   else
876                   {
877                     pastroppetit=Standard_False;
878                   }
879                 }
880               }
881               while(pastroppetit);
882             }
883
884             break;
885           }
886         case IntWalk_OK:
887         case IntWalk_ArretSurPoint://006
888           {
889             //=======================================================
890             //== Stop Test t   :  Frame on Param(.)     ==
891             //=======================================================
892             //xft arrive here
893             Arrive = TestArret(DejaReparti,Param,ChoixIso); 
894             // JMB 30th December 1999. 
895             // Some statement below should not be put in comment because they are useful.
896             // See grid CTO 909 A1 which infinitely loops 
897             if(Arrive==Standard_False && Status==IntWalk_ArretSurPoint)
898             {
899               Arrive=Standard_True;
900 #ifdef DEB
901               cout << "Compile with option DEB : if problems with intersection : ";
902               cout << "IntWalk_PWalking_1.gxx (lbr le 1erdec98)"<<endl;
903 #endif
904             }
905
906             if(Arrive)
907             {
908               NbPasOKConseq = -10;
909             }
910
911             if(!Arrive)//005
912             {
913               //=====================================================
914               //== Param(.) is in the limits                       ==
915               //==  and does not end a closed  line                ==
916               //=====================================================
917               //== Check on the current point of myInters
918               Standard_Boolean pointisvalid = Standard_False;
919               {
920                 Standard_Real u1,v1,u2,v2; 
921                 myIntersectionOn2S.Point().Parameters(u1,v1,u2,v2);
922
923                 //
924                 if(u1 <= UM1  && u2 <= UM2 && v1 <= VM1 && 
925                   v2 <= VM2  && u1 >= Um1 && u2 >= Um2 &&
926                   v1 >= Vm1  && v2 >= Vm2)
927                 {
928                   pointisvalid=Standard_True;
929                 }
930               }
931
932               //
933               if(pointisvalid)
934               {
935                 previousPoint = myIntersectionOn2S.Point();
936                 previoustg = myIntersectionOn2S.IsTangent();
937
938                 if(!previoustg)
939                 {
940                   previousd  = myIntersectionOn2S.Direction();
941                   previousd1 = myIntersectionOn2S.DirectionOnS1();
942                   previousd2 = myIntersectionOn2S.DirectionOnS2();
943                 }
944                 //=====================================================
945                 //== Check on the previous Point
946                 {
947                   Standard_Real u1,v1,u2,v2;
948                   previousPoint.Parameters(u1,v1,u2,v2); 
949                   if( u1 <= UM1  && u2 <= UM2 && v1 <= VM1 &&
950                     v2 <= VM2  && u1 >= Um1 && u2 >= Um2 &&
951                     v1 >= Vm1  && v2 >= Vm2)
952                   {
953                     pl = previousPoint.Value();
954                     if(bTestFirstPoint)
955                     {
956                       if(pf.SquareDistance(pl) < aSQDistMax)
957                       {
958                         IncKey++;
959                         if(IncKey == 5000)
960                           return;
961                         else
962                           continue;
963                       }
964                       else
965                       {
966                         bTestFirstPoint = Standard_False;
967                       }
968                     }
969                     //
970                     AddAPoint(line,previousPoint);
971                     RejectIndex++;
972
973                     if(RejectIndex >= RejectIndexMAX)
974                     {
975                       break;
976                     }
977
978                     //
979                     LevelOfIterWithoutAppend = 0;
980                   }
981                 }
982               }//pointisvalid
983               //====================================================
984
985               if(Status == IntWalk_ArretSurPoint)
986               {
987                 RepartirOuDiviser(DejaReparti,ChoixIso,Arrive);
988               }
989               else
990               {
991                 if (line->NbPoints() == 2)
992                 {
993                   pasSav[0] = pasuv[0];
994                   pasSav[1] = pasuv[1];
995                   pasSav[2] = pasuv[2];
996                   pasSav[3] = pasuv[3];
997                 }
998               }
999             }//005 if(!Arrive)
1000             else  //004
1001             {
1002               if(close)
1003               {
1004                 //================= la ligne est fermee ===============
1005                 AddAPoint(line,line->Value(1)); //ligne fermee
1006                 LevelOfIterWithoutAppend=0;
1007               }
1008               else    //$$$
1009               {
1010                 //====================================================
1011                 //== Param was not in the limits (was reframed)
1012                 //====================================================
1013                 Standard_Boolean bPrevNotTangent = !previoustg || !myIntersectionOn2S.IsTangent();
1014
1015                 IntImp_ConstIsoparametric SauvChoixIso = ChoixIso;
1016                 ChoixIso = myIntersectionOn2S.Perform(Param,Rsnld,ChoixIso);
1017                 //
1018                 if(!myIntersectionOn2S.IsEmpty()) //002
1019                 {
1020                   // mutially outpasses in the square or intersection in corner
1021
1022                   if(TestArret(Standard_True,Param,ChoixIso))
1023                   {
1024                     NbPasOKConseq = -10;
1025                     ChoixIso = myIntersectionOn2S.Perform(Param,Rsnld,ChoixIso);
1026
1027                     if(!myIntersectionOn2S.IsEmpty())
1028                     {
1029                       previousPoint = myIntersectionOn2S.Point();
1030                       previoustg = myIntersectionOn2S.IsTangent();
1031
1032                       if (!previoustg)
1033                       {
1034                         previousd  = myIntersectionOn2S.Direction();
1035                         previousd1 = myIntersectionOn2S.DirectionOnS1();
1036                         previousd2 = myIntersectionOn2S.DirectionOnS2();
1037                       }
1038
1039                       pl = previousPoint.Value();
1040
1041                       if(bTestFirstPoint)
1042                       {
1043                         if(pf.SquareDistance(pl) < aSQDistMax)
1044                         {
1045                           IncKey++;
1046                           if(IncKey == 5000)
1047                             return;
1048                           else
1049                             continue;
1050                         }
1051                         else
1052                         {
1053                           bTestFirstPoint = Standard_False;
1054                         }
1055                       }
1056                       //
1057                       AddAPoint(line,previousPoint);
1058                       RejectIndex++;
1059
1060                       if(RejectIndex >= RejectIndexMAX)
1061                       {
1062                         break;
1063                       }
1064
1065                       //
1066                       LevelOfIterWithoutAppend=0;
1067                       RepartirOuDiviser(DejaReparti,ChoixIso,Arrive);
1068                     }
1069                     else
1070                     {
1071                       //fail framing divides the step
1072                       Arrive = Standard_False;
1073                       RepartirOuDiviser(DejaReparti,ChoixIso,Arrive);
1074                       NoTestDeflection = Standard_True;
1075                       ChoixIso = SauvChoixIso;
1076                     }
1077                   }//if(TestArret())
1078                   else
1079                   {
1080                     // save the last point
1081                     // to revert to it if the current point is out of bounds
1082
1083                     IntSurf_PntOn2S previousPointSave = previousPoint;
1084                     Standard_Boolean previoustgSave   = previoustg;
1085                     gp_Dir previousdSave              = previousd;
1086                     gp_Dir2d previousd1Save           = previousd1;
1087                     gp_Dir2d previousd2Save           = previousd2;
1088
1089                     previousPoint = myIntersectionOn2S.Point();
1090                     previoustg = myIntersectionOn2S.IsTangent();
1091                     Arrive = Standard_False;
1092
1093                     if(!previoustg)
1094                     {
1095                       previousd  = myIntersectionOn2S.Direction();
1096                       previousd1 = myIntersectionOn2S.DirectionOnS1();
1097                       previousd2 = myIntersectionOn2S.DirectionOnS2();
1098                     }
1099
1100                     //========================================
1101                     //== Check on PreviousPoint @@
1102
1103                     {
1104                       Standard_Real u1,v1,u2,v2;
1105                       previousPoint.Parameters(u1,v1,u2,v2);
1106
1107                       //To save initial 2d points
1108                       gp_Pnt2d ParamPntOnS1(Param(1), Param(2));
1109                       gp_Pnt2d ParamPntOnS2(Param(3), Param(4));
1110
1111                       ///////////////////////////
1112                       Param(1) = u1;
1113                       Param(2) = v1;
1114                       Param(3) = u2;
1115                       Param(4) = v2;
1116                       //
1117
1118                       //xf
1119                       Standard_Boolean bFlag1, bFlag2;
1120                       Standard_Real aTol2D=1.e-11;
1121                       //
1122                       bFlag1=u1 >= Um1-aTol2D && v1 >= Vm1-aTol2D && u1 <= UM1+aTol2D && v1 <= VM1+aTol2D;
1123                       bFlag2=u2 >= Um2-aTol2D && v2 >= Vm2-aTol2D && u2 <= UM2+aTol2D && v2 <= VM2+aTol2D;
1124                       if (bFlag1 && bFlag2)
1125                       {
1126                         /*
1127                         if(u1 <= UM1  && u2 <= UM2 && v1 <= VM1 &&
1128                         v2 <= VM2  && u1 >= Um1 && u2 >= Um2 &&
1129                         v1 >= Vm1  && v2 >= Vm2)  {
1130                         */                      
1131                         //xt
1132                         pl = previousPoint.Value();
1133
1134                         if(bTestFirstPoint)
1135                         {
1136                           if(pf.SquareDistance(pl) < aSQDistMax)
1137                           {
1138                             IncKey++;
1139
1140                             if(IncKey == 5000)
1141                               return;
1142                             else
1143                               continue;
1144                           }
1145                           else
1146                           {
1147                             bTestFirstPoint = Standard_False;
1148                           }
1149                         }
1150
1151                         //To avoid walking around the same point
1152                         //in the tangent zone near a border
1153
1154                         if (previoustg)
1155                         {
1156                           Standard_Real prevU1, prevV1, prevU2, prevV2;
1157                           previousPointSave.Parameters(prevU1, prevV1, prevU2, prevV2);
1158                           gp_Pnt2d prevPntOnS1(prevU1, prevV1), prevPntOnS2(prevU2, prevV2);
1159                           gp_Pnt2d curPntOnS1(u1, v1), curPntOnS2(u2, v2);
1160                           gp_Vec2d PrevToParamOnS1(prevPntOnS1, ParamPntOnS1);
1161                           gp_Vec2d PrevToCurOnS1(prevPntOnS1, curPntOnS1);
1162                           gp_Vec2d PrevToParamOnS2(prevPntOnS2, ParamPntOnS2);
1163                           gp_Vec2d PrevToCurOnS2(prevPntOnS2, curPntOnS2);
1164                           Standard_Real MaxAngle = 3*M_PI/4;
1165
1166                           if (Abs(PrevToParamOnS1.Angle(PrevToCurOnS1)) > MaxAngle &&
1167                             Abs(PrevToParamOnS2.Angle(PrevToCurOnS2)) > MaxAngle)
1168                           {
1169                             Arrive = Standard_True;
1170                             break;
1171                           }
1172                         }
1173
1174                         ////////////////////////////////////////
1175                         AddAPoint(line,previousPoint);
1176                         RejectIndex++;
1177
1178                         if(RejectIndex >= RejectIndexMAX)
1179                         {
1180                           break;
1181                         }
1182
1183                         //
1184
1185                         LevelOfIterWithoutAppend=0;
1186                         Arrive = Standard_True;
1187                       }
1188                       else
1189                       {
1190                         // revert to the last correctly calculated point
1191                         previousPoint = previousPointSave;
1192                         previoustg    = previoustgSave;
1193                         previousd     = previousdSave;
1194                         previousd1    = previousd1Save;
1195                         previousd2    = previousd2Save;
1196                       }
1197                     }
1198
1199                     //
1200                     Standard_Boolean wasExtended = Standard_False;
1201
1202                     if(Arrive && myIntersectionOn2S.IsTangent() && bPrevNotTangent)
1203                     {
1204                       if(ExtendLineInCommonZone(SauvChoixIso, DejaReparti))
1205                       {
1206                         wasExtended = Standard_True;
1207                         Arrive = Standard_False;
1208                         ChoixIso = SauvChoixIso;
1209                       }
1210                     }
1211
1212                     RepartirOuDiviser(DejaReparti,ChoixIso,Arrive);
1213
1214                     if(Arrive && 
1215                       myIntersectionOn2S.IsDone() && !myIntersectionOn2S.IsEmpty() &&
1216                       myIntersectionOn2S.IsTangent() && bPrevNotTangent &&
1217                       !wasExtended)
1218                     {
1219                       if(ExtendLineInCommonZone(SauvChoixIso, DejaReparti))
1220                       {
1221                         wasExtended = Standard_True;
1222                         Arrive = Standard_False;
1223                         ChoixIso = SauvChoixIso;
1224                       }
1225                     }
1226                   }//else !TestArret() $
1227                 }//$$ end successful framing on border (!myIntersectionOn2S.IsEmpty())
1228                 else
1229                 {
1230                   //echec framing on border; division of step 
1231                   Arrive = Standard_False;
1232                   NoTestDeflection = Standard_True;
1233                   RepartirOuDiviser(DejaReparti,ChoixIso,Arrive);
1234                 }
1235               }//$$$ end framing on border (!close)
1236             }//004 fin TestArret return Arrive = True
1237           } // 006case IntWalk_ArretSurPoint:  end Processing Status = OK  or ArretSurPoint 
1238         } //007  switch(Status) 
1239       } //008 end processing point  (TEST DEFLECTION)
1240     } //009 end processing line (else if myIntersectionOn2S.IsDone())
1241   }  //010 end if first departure point allows marching  while (!Arrive)
1242
1243   done = Standard_True;
1244 }
1245 // ===========================================================================================================
1246 // function: ExtendLineInCommonZone
1247 // purpose:  Extends already computed line inside tangent zone in the direction given by theChoixIso.
1248 //           Returns Standard_True if the line was extended through tangent zone and the last computed point 
1249 //           is outside the tangent zone (but it is not put into the line). Otherwise returns Standard_False.
1250 // ===========================================================================================================
1251 Standard_Boolean IntWalk_PWalking::ExtendLineInCommonZone(const IntImp_ConstIsoparametric theChoixIso,
1252                                                           const Standard_Boolean          theDirectionFlag) 
1253 {
1254   Standard_Boolean bOutOfTangentZone = Standard_False;
1255   Standard_Boolean bStop = !myIntersectionOn2S.IsTangent();
1256   Standard_Integer dIncKey = 1;
1257   TColStd_Array1OfReal Param(1,4);
1258   IntWalk_StatusDeflection Status = IntWalk_OK;
1259   Standard_Integer nbIterWithoutAppend = 0;
1260   Standard_Integer nbEqualPoints = 0;
1261   Standard_Integer parit = 0;
1262   Standard_Integer uvit = 0;
1263   IntSurf_SequenceOfPntOn2S aSeqOfNewPoint;
1264
1265   while (!bStop) {
1266     nbIterWithoutAppend++;
1267
1268     if((nbIterWithoutAppend > 20) || (nbEqualPoints > 20)) {
1269 #ifdef DEB
1270       cout<<"Compile with option DEB:";
1271       cout<<"Infinite loop has detected. Stop iterations (IntWalk_PWalking_1.gxx)" << endl;
1272 #endif
1273       bStop = Standard_True;
1274       break;
1275     }
1276     Standard_Real f = 0.;
1277
1278     switch (theChoixIso)
1279     { 
1280     case IntImp_UIsoparametricOnCaro1: f = Abs(previousd1.X()); break;
1281     case IntImp_VIsoparametricOnCaro1: f = Abs(previousd1.Y()); break;
1282     case IntImp_UIsoparametricOnCaro2: f = Abs(previousd2.X()); break;
1283     case IntImp_VIsoparametricOnCaro2: f = Abs(previousd2.Y()); break;
1284     }
1285
1286     if(f<0.1) f=0.1;
1287
1288     previousPoint.Parameters(Param(1),Param(2),Param(3),Param(4));
1289
1290     Standard_Real dP1 = sensCheminement * pasuv[0] * previousd1.X() /f;
1291     Standard_Real dP2 = sensCheminement * pasuv[1] * previousd1.Y() /f;
1292     Standard_Real dP3 = sensCheminement * pasuv[2] * previousd2.X() /f; 
1293     Standard_Real dP4 = sensCheminement * pasuv[3] * previousd2.Y() /f;
1294
1295     if(theChoixIso == IntImp_UIsoparametricOnCaro1 && Abs(dP1) < 1.e-7) dP1 *= (5. * (Standard_Real)dIncKey);
1296     if(theChoixIso == IntImp_VIsoparametricOnCaro1 && Abs(dP2) < 1.e-7) dP2 *= (5. * (Standard_Real)dIncKey);
1297     if(theChoixIso == IntImp_UIsoparametricOnCaro2 && Abs(dP3) < 1.e-7) dP3 *= (5. * (Standard_Real)dIncKey);
1298     if(theChoixIso == IntImp_VIsoparametricOnCaro2 && Abs(dP4) < 1.e-7) dP4 *= (5. * (Standard_Real)dIncKey);
1299
1300     Param(1) += dP1;
1301     Param(2) += dP2;
1302     Param(3) += dP3; 
1303     Param(4) += dP4;
1304     Standard_Real SvParam[4];
1305     IntImp_ConstIsoparametric ChoixIso = theChoixIso;
1306
1307     for(parit = 0; parit < 4; parit++) {
1308       SvParam[parit] = Param(parit+1);
1309     }
1310     math_FunctionSetRoot  Rsnld(myIntersectionOn2S.Function());
1311     ChoixIso = myIntersectionOn2S.Perform(Param,Rsnld, theChoixIso);
1312
1313     if (!myIntersectionOn2S.IsDone()) {
1314       return bOutOfTangentZone;
1315     }
1316     else {
1317       if (myIntersectionOn2S.IsEmpty()) {
1318         return bOutOfTangentZone;
1319       }
1320
1321       Status = TestDeflection();
1322
1323       if(Status == IntWalk_OK) {
1324
1325         for(uvit = 0; uvit < 4; uvit++) {
1326           if(pasuv[uvit] < pasInit[uvit]) {
1327             pasuv[uvit] = pasInit[uvit];
1328           }
1329         }
1330       }
1331
1332       switch(Status) {
1333       case  IntWalk_ArretSurPointPrecedent:
1334         {
1335           bStop = Standard_True;
1336           bOutOfTangentZone = !myIntersectionOn2S.IsTangent();
1337           break;
1338         }
1339       case IntWalk_PasTropGrand:
1340         {
1341           for(parit = 0; parit < 4; parit++) {
1342             Param(parit+1) = SvParam[parit];
1343           }
1344           Standard_Boolean bDecrease = Standard_False;
1345
1346           for(uvit = 0; uvit < 4; uvit++) {
1347             if(pasSav[uvit] < pasInit[uvit]) { 
1348               pasInit[uvit] -= (pasInit[uvit] - pasSav[uvit]) * 0.1;
1349               bDecrease = Standard_True;
1350             }
1351           }
1352
1353           if(bDecrease) nbIterWithoutAppend--;
1354           break;
1355         }
1356       case IntWalk_PointConfondu:
1357         {
1358           for(uvit = 0; uvit < 4; uvit++) {
1359             if(pasuv[uvit] < pasInit[uvit]) {
1360               pasuv[uvit] += (pasInit[uvit] - pasuv[uvit]) * 0.1;
1361             }
1362           }
1363           break;
1364         }
1365       case IntWalk_OK:
1366       case IntWalk_ArretSurPoint:
1367         {
1368           //
1369           bStop = TestArret(theDirectionFlag, Param, ChoixIso);
1370           //
1371
1372           //
1373           if(!bStop) {
1374             Standard_Real u11,v11,u12,v12; 
1375             myIntersectionOn2S.Point().Parameters(u11,v11,u12,v12); 
1376             Standard_Real u21,v21,u22,v22;
1377             previousPoint.Parameters(u21,v21,u22,v22); 
1378
1379             if(((fabs(u11-u21) < ResoU1) && (fabs(v11-v21) < ResoV1)) ||
1380               ((fabs(u12-u22) < ResoU2) && (fabs(v12-v22) < ResoV2))) {
1381                 nbEqualPoints++;
1382             }
1383             else {
1384               nbEqualPoints = 0;
1385             }
1386           }
1387           //
1388
1389           bStop = bStop || !myIntersectionOn2S.IsTangent();
1390           bOutOfTangentZone = !myIntersectionOn2S.IsTangent();
1391
1392           if(!bStop) {
1393             Standard_Boolean pointisvalid = Standard_False;
1394             Standard_Real u1,v1,u2,v2; 
1395             myIntersectionOn2S.Point().Parameters(u1,v1,u2,v2); 
1396
1397             if(u1 <= UM1  && u2 <= UM2 && v1 <= VM1 && 
1398               v2 <= VM2  && u1 >= Um1 && u2 >= Um2 &&
1399               v1 >= Vm1  && v2 >= Vm2) 
1400               pointisvalid = Standard_True;
1401
1402             if(pointisvalid) {
1403               previousPoint = myIntersectionOn2S.Point();
1404               previoustg = myIntersectionOn2S.IsTangent();
1405
1406               if(!previoustg) {
1407                 previousd  = myIntersectionOn2S.Direction();
1408                 previousd1 = myIntersectionOn2S.DirectionOnS1();
1409                 previousd2 = myIntersectionOn2S.DirectionOnS2();
1410               }
1411               Standard_Boolean bAddPoint = Standard_True;
1412
1413               if(line->NbPoints() >= 1) {
1414                 gp_Pnt pf = line->Value(1).Value();
1415                 gp_Pnt pl = previousPoint.Value(); 
1416
1417                 if(pf.Distance(pl) < Precision::Confusion()) { 
1418                   dIncKey++; 
1419                   if(dIncKey == 5000) return bOutOfTangentZone; 
1420                   else bAddPoint = Standard_False;
1421                 }
1422               }
1423
1424               if(bAddPoint) {
1425                 aSeqOfNewPoint.Append(previousPoint);
1426                 nbIterWithoutAppend = 0;
1427               }
1428             }
1429
1430             if (line->NbPoints() == 2) {
1431               for(uvit = 0; uvit < 4; uvit++) {
1432                 pasSav[uvit] = pasuv[uvit]; 
1433               }
1434             }
1435
1436             if ( !pointisvalid ) {
1437               // decrease step if out of bounds
1438               // otherwise the same calculations will be 
1439               // repeated several times
1440               if ( ( u1 > UM1 ) || ( u1 < Um1 ) )
1441                 pasuv[0] *= 0.5;
1442
1443               if ( ( v1 > VM1 ) || ( v1 < Vm1 ) ) 
1444                 pasuv[1] *= 0.5;
1445
1446               if ( ( u2 > UM2 ) || ( u2 < Um2 ) )
1447                 pasuv[2] *= 0.5;
1448
1449               if ( ( v2 > VM2 ) || ( v2 < Vm2 ) )
1450                 pasuv[3] *= 0.5;
1451             }
1452           } // end if(!bStop)
1453           else { //if(bStop)
1454             if(close && (line->NbPoints() >= 1)) { 
1455
1456               if(!bOutOfTangentZone) {
1457                 aSeqOfNewPoint.Append(line->Value(1)); // line end
1458               }
1459               nbIterWithoutAppend = 0;
1460             }
1461             else {
1462               ChoixIso = myIntersectionOn2S.Perform(Param, Rsnld, theChoixIso);
1463
1464               if(myIntersectionOn2S.IsEmpty()) { 
1465                 bStop = !myIntersectionOn2S.IsTangent();
1466                 bOutOfTangentZone = !myIntersectionOn2S.IsTangent();
1467               }
1468               else {
1469                 Standard_Boolean bAddPoint = Standard_True;
1470                 Standard_Boolean pointisvalid = Standard_False;
1471
1472                 previousPoint = myIntersectionOn2S.Point();
1473                 Standard_Real u1,v1,u2,v2; 
1474                 previousPoint.Parameters(u1,v1,u2,v2); 
1475
1476                 if(u1 <= UM1  && u2 <= UM2 && v1 <= VM1 && 
1477                   v2 <= VM2  && u1 >= Um1 && u2 >= Um2 &&
1478                   v1 >= Vm1  && v2 >= Vm2) 
1479                   pointisvalid = Standard_True;
1480
1481                 if(pointisvalid) {
1482
1483                   if(line->NbPoints() >= 1) {
1484                     gp_Pnt pf = line->Value(1).Value();
1485                     gp_Pnt pl = previousPoint.Value(); 
1486
1487                     if(pf.Distance(pl) < Precision::Confusion()) { 
1488                       dIncKey++; 
1489                       if(dIncKey == 5000) return bOutOfTangentZone; 
1490                       else bAddPoint = Standard_False;
1491                     }
1492                   }
1493
1494                   if(bAddPoint && !bOutOfTangentZone) {
1495                     aSeqOfNewPoint.Append(previousPoint);
1496                     nbIterWithoutAppend = 0;
1497                   }
1498                 }
1499               }
1500             }
1501           }
1502           break;
1503         }
1504       default:
1505         {
1506           break;
1507         }
1508       }
1509     }
1510   }
1511   Standard_Boolean bExtendLine = Standard_False;
1512   Standard_Real u1 = 0., v1 = 0., u2 = 0., v2 = 0.; 
1513
1514   Standard_Integer pit = 0;
1515
1516   for(pit = 0; !bExtendLine && (pit < 2); pit++) {
1517     if(pit == 0)
1518       previousPoint.Parameters(u1,v1,u2,v2); 
1519     else {
1520       if(aSeqOfNewPoint.Length() > 0)
1521         aSeqOfNewPoint.Value(aSeqOfNewPoint.Length()).Parameters(u1,v1,u2,v2); 
1522       else
1523         break;
1524     }
1525
1526     if(((u1 - Um1) < ResoU1) ||
1527       ((UM1 - u1) < ResoU1) ||
1528       ((u2 - Um2) < ResoU2) ||
1529       ((UM2 - u2) < ResoU2) ||
1530       ((v1 - Vm1) < ResoV1) ||
1531       ((VM1 - v1) < ResoV1) ||
1532       ((v2 - Vm2) < ResoV2) ||
1533       ((VM2 - v2) < ResoV2))
1534       bExtendLine = Standard_True;
1535   }
1536
1537   if(!bExtendLine) {
1538     //    if(Status == IntWalk_OK || Status == IntWalk_ArretSurPoint) {
1539     if(Status == IntWalk_OK) {
1540       bExtendLine = Standard_True;
1541
1542       if(aSeqOfNewPoint.Length() > 1) {
1543         TColStd_Array1OfReal FirstParams(0, 3), LastParams(0, 3), Resolutions(0, 3);
1544         Resolutions(0) = ResoU1; Resolutions(1) = ResoV1; Resolutions(2) = ResoU2; Resolutions(3) = ResoV2;
1545
1546         aSeqOfNewPoint(1).Parameters(FirstParams.ChangeValue(0), FirstParams.ChangeValue(1),
1547           FirstParams.ChangeValue(2), FirstParams.ChangeValue(3));
1548         aSeqOfNewPoint(aSeqOfNewPoint.Length()).Parameters(LastParams.ChangeValue(0), 
1549           LastParams.ChangeValue(1),
1550           LastParams.ChangeValue(2), 
1551           LastParams.ChangeValue(3)); 
1552         Standard_Integer indexofiso = 0;
1553
1554         if(theChoixIso == IntImp_UIsoparametricOnCaro1) indexofiso = 0;
1555         if(theChoixIso == IntImp_VIsoparametricOnCaro1) indexofiso = 1;
1556         if(theChoixIso == IntImp_UIsoparametricOnCaro2) indexofiso = 2;
1557         if(theChoixIso == IntImp_VIsoparametricOnCaro2) indexofiso = 3;
1558
1559         Standard_Integer afirstindex = (indexofiso < 2) ? 0 : 2;
1560         gp_Vec2d aTangentZoneDir(gp_Pnt2d(FirstParams.Value(afirstindex), FirstParams.Value(afirstindex + 1)),
1561           gp_Pnt2d(LastParams.Value(afirstindex), LastParams.Value(afirstindex + 1)));
1562
1563         gp_Dir2d anIsoDir(0, 1);
1564
1565         if((indexofiso == 1) || (indexofiso == 3))
1566           anIsoDir = gp_Dir2d(1, 0);
1567
1568         if(aTangentZoneDir.SquareMagnitude() > gp::Resolution()) {
1569           Standard_Real piquota = M_PI*0.25;
1570
1571           if(fabs(aTangentZoneDir.Angle(anIsoDir)) > piquota) {
1572             Standard_Integer ii = 1, nextii = 2;
1573             gp_Vec2d d1(0, 0);
1574             Standard_Real asqresol = gp::Resolution();
1575             asqresol *= asqresol;
1576
1577             do {
1578               aSeqOfNewPoint(ii).Parameters(FirstParams.ChangeValue(0), FirstParams.ChangeValue(1),
1579                 FirstParams.ChangeValue(2), FirstParams.ChangeValue(3));
1580               aSeqOfNewPoint(ii + 1).Parameters(LastParams.ChangeValue(0), LastParams.ChangeValue(1),
1581                 LastParams.ChangeValue(2), LastParams.ChangeValue(3));
1582               d1 = gp_Vec2d(gp_Pnt2d(FirstParams.Value(afirstindex),
1583                 FirstParams.Value(afirstindex + 1)),
1584                 gp_Pnt2d(LastParams.Value(afirstindex),
1585                 LastParams.Value(afirstindex + 1)));
1586               ii++;
1587             }
1588             while((d1.SquareMagnitude() < asqresol) &&
1589               (ii < aSeqOfNewPoint.Length()));
1590
1591             nextii = ii;
1592
1593             while(nextii < aSeqOfNewPoint.Length()) {
1594
1595               gp_Vec2d nextd1(0, 0);
1596               Standard_Integer jj = nextii;
1597
1598               do {
1599                 aSeqOfNewPoint(jj).Parameters(FirstParams.ChangeValue(0), FirstParams.ChangeValue(1),
1600                   FirstParams.ChangeValue(2), FirstParams.ChangeValue(3));
1601                 aSeqOfNewPoint(jj + 1).Parameters(LastParams.ChangeValue(0), LastParams.ChangeValue(1),
1602                   LastParams.ChangeValue(2), LastParams.ChangeValue(3));
1603                 nextd1 = gp_Vec2d(gp_Pnt2d(FirstParams.Value(afirstindex),
1604                   FirstParams.Value(afirstindex + 1)),
1605                   gp_Pnt2d(LastParams.Value(afirstindex),
1606                   LastParams.Value(afirstindex + 1)));
1607                 jj++;
1608
1609               }
1610               while((nextd1.SquareMagnitude() < asqresol) &&
1611                 (jj < aSeqOfNewPoint.Length()));
1612               nextii = jj;
1613
1614               if(fabs(d1.Angle(nextd1)) > piquota) {
1615                 bExtendLine = Standard_False;
1616                 break;
1617               }
1618               d1 = nextd1;
1619             }
1620           }
1621           // end if(fabs(aTangentZoneDir.Angle(anIsoDir)
1622         }
1623       }
1624     }
1625   }
1626
1627   if(!bExtendLine) {
1628     return Standard_False;
1629   }
1630   Standard_Integer i = 0;
1631
1632   for(i = 1; i <= aSeqOfNewPoint.Length(); i++) {
1633     AddAPoint(line, aSeqOfNewPoint.Value(i));
1634   }
1635
1636   return bOutOfTangentZone;
1637 }
1638
1639 Standard_Boolean IntWalk_PWalking::DistanceMinimizeByGradient(const Handle(Adaptor3d_HSurface)& theASurf1,
1640                                                               const Handle(Adaptor3d_HSurface)& theASurf2,
1641                                                               Standard_Real& theU1,
1642                                                               Standard_Real& theV1,
1643                                                               Standard_Real& theU2,
1644                                                               Standard_Real& theV2,
1645                                                               const Standard_Real theStep0U1V1,
1646                                                               const Standard_Real theStep0U2V2)
1647 {
1648   const Standard_Integer aNbIterMAX = 60;
1649   const Standard_Real aTol = 1.0e-14;
1650   Handle(Geom_Surface) aS1, aS2;
1651
1652   switch(theASurf1->GetType())
1653   {
1654   case GeomAbs_BezierSurface:
1655     aS1 = theASurf1->Surface().Bezier();
1656     break;
1657   case GeomAbs_BSplineSurface:
1658     aS1 = theASurf1->Surface().BSpline();
1659     break;
1660   default:
1661     return Standard_True;
1662   }
1663
1664   switch(theASurf2->GetType())
1665   {
1666   case GeomAbs_BezierSurface:
1667     aS2 = theASurf2->Surface().Bezier();
1668     break;
1669   case GeomAbs_BSplineSurface:
1670     aS2 = theASurf2->Surface().BSpline();
1671     break;
1672   default:
1673     return Standard_True;
1674   }
1675
1676   Standard_Boolean aStatus = Standard_False;
1677
1678   gp_Pnt aP1, aP2;
1679   gp_Vec aD1u, aD1v, aD2U, aD2V;
1680
1681   aS1->D1(theU1, theV1, aP1, aD1u, aD1v);
1682   aS2->D1(theU2, theV2, aP2, aD2U, aD2V);
1683
1684   Standard_Real aSQDistPrev = aP1.SquareDistance(aP2);
1685
1686   gp_Vec aP12(aP1, aP2);
1687
1688   Standard_Real aGradFu(-aP12.Dot(aD1u));
1689   Standard_Real aGradFv(-aP12.Dot(aD1v));
1690   Standard_Real aGradFU( aP12.Dot(aD2U));
1691   Standard_Real aGradFV( aP12.Dot(aD2V));
1692
1693   Standard_Real aSTEPuv = theStep0U1V1, aStepUV = theStep0U2V2;
1694
1695   Standard_Boolean flRepeat = Standard_True;
1696   Standard_Integer aNbIter = aNbIterMAX;
1697
1698   while(flRepeat)
1699   {
1700     Standard_Real anAdd = aGradFu*aSTEPuv;
1701     Standard_Real aPARu = (anAdd >= 0.0)? (theU1 - Max(anAdd, Epsilon(theU1))) : (theU1 + Max(-anAdd, Epsilon(theU1)));
1702     anAdd = aGradFv*aSTEPuv;
1703     Standard_Real aPARv = (anAdd >= 0.0)? (theV1 - Max(anAdd, Epsilon(theV1))) : (theV1 + Max(-anAdd, Epsilon(theV1)));
1704     anAdd = aGradFU*aStepUV;
1705     Standard_Real aParU = (anAdd >= 0.0)? (theU2 - Max(anAdd, Epsilon(theU2))) : (theU2 + Max(-anAdd, Epsilon(theU2)));
1706     anAdd = aGradFV*aStepUV;
1707     Standard_Real aParV = (anAdd >= 0.0)? (theV2 - Max(anAdd, Epsilon(theV2))) : (theV2 + Max(-anAdd, Epsilon(theV2)));
1708
1709     gp_Pnt aPt1, aPt2;
1710
1711     aS1->D1(aPARu, aPARv, aPt1, aD1u, aD1v);
1712     aS2->D1(aParU, aParV, aPt2, aD2U, aD2V);
1713
1714     Standard_Real aSQDist = aPt1.SquareDistance(aPt2);
1715
1716     if(aSQDist < aSQDistPrev)
1717     {
1718       aSQDistPrev = aSQDist;
1719       theU1 = aPARu;
1720       theV1 = aPARv;
1721       theU2 = aParU;
1722       theV2 = aParV;
1723
1724       aStatus = aSQDistPrev < aTol;
1725       aSTEPuv *= 1.2;
1726       aStepUV *= 1.2;
1727     }
1728     else
1729     {
1730       if(--aNbIter < 0)
1731       {
1732         flRepeat = Standard_False;
1733       }
1734       else
1735       {
1736         aS1->D1(theU1, theV1, aPt1, aD1u, aD1v);
1737         aS2->D1(theU2, theV2, aPt2, aD2U, aD2V);
1738
1739         gp_Vec aP12(aPt1, aPt2);
1740         aGradFu = -aP12.Dot(aD1u);
1741         aGradFv = -aP12.Dot(aD1v);
1742         aGradFU = aP12.Dot(aD2U);
1743         aGradFV = aP12.Dot(aD2V);
1744         aSTEPuv = theStep0U1V1;
1745         aStepUV = theStep0U2V2;
1746       }
1747     }
1748   }
1749
1750   return aStatus;
1751 }
1752
1753 Standard_Boolean IntWalk_PWalking::DistanceMinimizeByExtrema( const Handle(Adaptor3d_HSurface)& theASurf, 
1754                                                               const gp_Pnt& theP0,
1755                                                               Standard_Real& theU0,
1756                                                               Standard_Real& theV0,
1757                                                               const Standard_Real theStep0U,
1758                                                               const Standard_Real theStep0V)
1759 {
1760   const Standard_Real aTol = 1.0e-14;
1761   gp_Pnt aPS;
1762   gp_Vec aD1Su, aD1Sv, aD2Su, aD2Sv, aD2SuvTemp;
1763   Standard_Real aSQDistPrev = RealLast();
1764   Standard_Real aU = theU0, aV = theV0;
1765   
1766   Standard_Integer aNbIter = 10;
1767   do
1768   {
1769     theASurf->D2(aU, aV, aPS, aD1Su, aD1Sv, aD2Su, aD2Sv, aD2SuvTemp);
1770     
1771     gp_Vec aVec(theP0, aPS);
1772     
1773     Standard_Real aSQDist = aVec.SquareMagnitude();
1774
1775     if(aSQDist >= aSQDistPrev)
1776       break;
1777
1778     aSQDistPrev = aSQDist;
1779     theU0 = aU;
1780     theV0 = aV;
1781     aNbIter--;
1782
1783     if(aSQDistPrev < aTol)
1784       break;
1785
1786     //Functions
1787     const Standard_Real aF1 = aD1Su.Dot(aVec), aF2 = aD1Sv.Dot(aVec);
1788
1789     //Derivatives
1790     const Standard_Real aDf1u = aD2Su.Dot(aVec) + aD1Su.Dot(aD1Su),
1791                         aDf1v = aD2Su.Dot(aD1Sv),
1792                         aDf2u = aDf1v,
1793                         aDf2v = aD2Sv.Dot(aVec) + aD1Sv.Dot(aD1Sv);
1794
1795     const Standard_Real aDet = aDf1u*aDf2v - aDf1v*aDf2u;
1796     aU -= theStep0U*(aDf2v*aF1 - aDf1v*aF2)/aDet;
1797     aV += theStep0V*(aDf2u*aF1 - aDf1u*aF2)/aDet;
1798   }
1799   while(aNbIter > 0);
1800
1801   return (aSQDistPrev < aTol);
1802 }
1803
1804 Standard_Boolean IntWalk_PWalking::SeekAdditionalPoints(const Handle(Adaptor3d_HSurface)& theASurf1,
1805                                                         const Handle(Adaptor3d_HSurface)& theASurf2,
1806                                                         const Standard_Integer theMinNbPoints)
1807 {
1808   const Standard_Real aTol = 1.0e-14;
1809   Standard_Integer aNbPoints = line->NbPoints();
1810   if(aNbPoints > theMinNbPoints)
1811     return Standard_True;
1812
1813   const Standard_Real aU1bFirst = theASurf1->FirstUParameter();
1814   const Standard_Real aU1bLast = theASurf1->LastUParameter();
1815   const Standard_Real aU2bFirst = theASurf2->FirstUParameter();
1816   const Standard_Real aU2bLast = theASurf2->LastUParameter();
1817   const Standard_Real aV1bFirst = theASurf1->FirstVParameter();
1818   const Standard_Real aV1bLast = theASurf1->LastVParameter();
1819   const Standard_Real aV2bFirst = theASurf2->FirstVParameter();
1820   const Standard_Real aV2bLast = theASurf2->LastVParameter();
1821
1822   
1823   Standard_Boolean isPrecise = Standard_False;
1824
1825   Standard_Real U1prec = 0.0, V1prec = 0.0, U2prec = 0.0, V2prec = 0.0;
1826
1827   Standard_Integer aNbPointsPrev = 0;
1828   while(aNbPoints < theMinNbPoints && (aNbPoints != aNbPointsPrev))
1829   {
1830     aNbPointsPrev = aNbPoints;
1831     for(Standard_Integer fp = 1, lp = 2; fp < aNbPoints; fp = lp + 1)
1832     {
1833       Standard_Real U1f, V1f, U2f, V2f; //first point in 1st and 2nd surafaces
1834       Standard_Real U1l, V1l, U2l, V2l; //last  point in 1st and 2nd surafaces
1835
1836       lp = fp+1;
1837       line->Value(fp).Parameters(U1f, V1f, U2f, V2f);
1838       line->Value(lp).Parameters(U1l, V1l, U2l, V2l);
1839
1840       U1prec = 0.5*(U1f+U1l);
1841       if(U1prec < aU1bFirst)
1842         U1prec = aU1bFirst;
1843       if(U1prec > aU1bLast)
1844         U1prec = aU1bLast;
1845
1846       V1prec = 0.5*(V1f+V1l);
1847       if(V1prec < aV1bFirst)
1848         V1prec = aV1bFirst;
1849       if(V1prec > aV1bLast)
1850         V1prec = aV1bLast;
1851
1852       U2prec = 0.5*(U2f+U2l);
1853       if(U2prec < aU2bFirst)
1854         U2prec = aU2bFirst;
1855       if(U2prec > aU2bLast)
1856         U2prec = aU2bLast;
1857
1858       V2prec = 0.5*(V2f+V2l);
1859       if(V2prec < aV2bFirst)
1860         V2prec = aV2bFirst;
1861       if(V2prec > aV2bLast)
1862         V2prec = aV2bLast;
1863
1864       Standard_Boolean aStatus = Standard_False;
1865       Standard_Integer aNbIter = 5;
1866       do
1867       {
1868         aStatus = DistanceMinimizeByGradient(theASurf1, theASurf2, U1prec, V1prec, U2prec, V2prec);
1869         if(aStatus)
1870         {
1871           break;
1872         }
1873
1874         aStatus = DistanceMinimizeByExtrema(theASurf1, theASurf2->Value(U2prec, V2prec), U1prec, V1prec);
1875         if(aStatus)
1876         {
1877           break;
1878         }
1879
1880         aStatus = DistanceMinimizeByExtrema(theASurf2, theASurf1->Value(U1prec, V1prec), U2prec, V2prec);
1881         if(aStatus)
1882         {
1883           break;
1884         }
1885       }
1886       while(!aStatus && (--aNbIter > 0));
1887
1888       if(aStatus)
1889       {
1890         gp_Pnt  aP1 = theASurf1->Value(U1prec, V1prec),
1891                 aP2 = theASurf2->Value(U2prec, V2prec);
1892         gp_Pnt aPInt(0.5*(aP1.XYZ() + aP2.XYZ()));
1893
1894         const Standard_Real aSQDist1 = aPInt.SquareDistance(aP1),
1895                             aSQDist2 = aPInt.SquareDistance(aP2);
1896
1897         if((aSQDist1 < aTol) && (aSQDist2 < aTol))
1898         {
1899           IntSurf_PntOn2S anIP;
1900           anIP.SetValue(aPInt, U1prec, V1prec, U2prec, V2prec);
1901           line->InsertBefore(lp, anIP);
1902           
1903           isPrecise = Standard_True;
1904
1905           if(++aNbPoints >= theMinNbPoints)
1906             break;
1907         }
1908         else
1909         {
1910           lp--;
1911         }
1912       }
1913     }
1914   }
1915
1916   return isPrecise;
1917 }