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