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