src/IntTools/IntTools_FaceFace.cxx

   1 // Created on: 2000-11-23
   2 // Created by: Michael KLOKOV
   3 // Copyright (c) 2000-2014 OPEN CASCADE SAS
   4 //
   5 // This file is part of Open CASCADE Technology software library.
   6 //
   7 // This library is free software; you can redistribute it and/or modify it under
   8 // the terms of the GNU Lesser General Public License version 2.1 as published
   9 // by the Free Software Foundation, with special exception defined in the file
  10 // OCCT_LGPL_EXCEPTION.txt. Consult the file LICENSE_LGPL_21.txt included in OCCT
  11 // distribution for complete text of the license and disclaimer of any warranty.
  12 //
  13 // Alternatively, this file may be used under the terms of Open CASCADE
  14 // commercial license or contractual agreement.
  15
  16 #include <IntTools_FaceFace.hxx>
  17
  18 #include <BRepAdaptor_Surface.hxx>
  19 #include <BRepTools.hxx>
  20 #include <BRep_Tool.hxx>
  21 #include <ElCLib.hxx>
  22 #include <ElSLib.hxx>
  23 #include <Geom2dAdaptor_Curve.hxx>
  24 #include <Geom2dInt_GInter.hxx>
  25 #include <Geom2d_Line.hxx>
  26 #include <Geom2d_TrimmedCurve.hxx>
  27 #include <GeomAPI_ProjectPointOnSurf.hxx>
  28 #include <GeomAdaptor_HSurface.hxx>
  29 #include <GeomInt_IntSS.hxx>
  30 #include <GeomInt_WLApprox.hxx>
  31 #include <GeomLib_Check2dBSplineCurve.hxx>
  32 #include <GeomLib_CheckBSplineCurve.hxx>
  33 #include <Geom_BSplineCurve.hxx>
  34 #include <Geom_Circle.hxx>
  35 #include <Geom_Ellipse.hxx>
  36 #include <Geom_Hyperbola.hxx>
  37 #include <Geom_Line.hxx>
  38 #include <Geom_OffsetSurface.hxx>
  39 #include <Geom_Parabola.hxx>
  40 #include <Geom_RectangularTrimmedSurface.hxx>
  41 #include <Geom_TrimmedCurve.hxx>
  42 #include <IntAna_QuadQuadGeo.hxx>
  43 #include <IntPatch_GLine.hxx>
  44 #include <IntPatch_RLine.hxx>
  45 #include <IntPatch_WLine.hxx>
  46 #include <IntRes2d_Domain.hxx>
  47 #include <IntSurf_Quadric.hxx>
  48 #include <IntTools_Context.hxx>
  49 #include <IntTools_Tools.hxx>
  50 #include <IntTools_TopolTool.hxx>
  51 #include <IntTools_WLineTool.hxx>
  52 #include <ProjLib_Plane.hxx>
  53 #include <TopExp_Explorer.hxx>
  54 #include <TopoDS.hxx>
  55 #include <TopoDS_Edge.hxx>
  56 #include <gp_Elips.hxx>
  57
  58 static
  59   void TolR3d(const Standard_Real aTolF1,
  60               const Standard_Real aTolF2,
  61               Standard_Real& myTolReached3d);
  62
  63 static
  64   void Parameters(const Handle(GeomAdaptor_HSurface)&,
  65                   const Handle(GeomAdaptor_HSurface)&,
  66                   const gp_Pnt&,
  67                   Standard_Real&,
  68                   Standard_Real&,
  69                   Standard_Real&,
  70                   Standard_Real&);
  71
  72 static
  73   void CorrectSurfaceBoundaries(const TopoDS_Face&  theFace,
  74                                 const Standard_Real theTolerance,
  75                                 Standard_Real&      theumin,
  76                                 Standard_Real&      theumax,
  77                                 Standard_Real&      thevmin,
  78                                 Standard_Real&      thevmax);
  79
  80 static
  81   Standard_Boolean ParameterOutOfBoundary(const Standard_Real       theParameter,
  82                                           const Handle(Geom_Curve)& theCurve,
  83                                           const TopoDS_Face&        theFace1,
  84                                           const TopoDS_Face&        theFace2,
  85                                           const Standard_Real       theOtherParameter,
  86                                           const Standard_Boolean    bIncreasePar,
  87                                           const Standard_Real       theTol,
  88                                           Standard_Real&            theNewParameter,
  89                                           const Handle(IntTools_Context)& );
  90
  91 static
  92   Standard_Boolean IsCurveValid(const Handle(Geom2d_Curve)& thePCurve);
  93
  94 static
  95   Standard_Boolean  ApproxWithPCurves(const gp_Cylinder& theCyl,
  96                                       const gp_Sphere& theSph);
  97
  98 static void  PerformPlanes(const Handle(GeomAdaptor_HSurface)& theS1,
  99                            const Handle(GeomAdaptor_HSurface)& theS2,
 100                            const Standard_Real TolF1,
 101                            const Standard_Real TolF2,
 102                            const Standard_Real TolAng,
 103                            const Standard_Real TolTang,
 104                            const Standard_Boolean theApprox1,
 105                            const Standard_Boolean theApprox2,
 106                            IntTools_SequenceOfCurves& theSeqOfCurve,
 107                            Standard_Boolean& theTangentFaces,
 108                            Standard_Real& TolReached3d);
 109
 110 static Standard_Boolean ClassifyLin2d(const Handle(GeomAdaptor_HSurface)& theS,
 111                                       const gp_Lin2d& theLin2d,
 112                                       const Standard_Real theTol,
 113                                       Standard_Real& theP1,
 114                                       Standard_Real& theP2);
 115 //
 116 static
 117   void ApproxParameters(const Handle(GeomAdaptor_HSurface)& aHS1,
 118                         const Handle(GeomAdaptor_HSurface)& aHS2,
 119                         Standard_Integer& iDegMin,
 120                         Standard_Integer& iNbIter,
 121                         Standard_Integer& iDegMax);
 122
 123 static
 124   void Tolerances(const Handle(GeomAdaptor_HSurface)& aHS1,
 125                   const Handle(GeomAdaptor_HSurface)& aHS2,
 126                   Standard_Real& aTolTang);
 127
 128 static
 129   Standard_Boolean SortTypes(const GeomAbs_SurfaceType aType1,
 130                              const GeomAbs_SurfaceType aType2);
 131 static
 132   Standard_Integer IndexType(const GeomAbs_SurfaceType aType);
 133
 134 //
 135 static
 136   Standard_Boolean CheckPCurve(const Handle(Geom2d_Curve)& aPC,
 137                                const TopoDS_Face& aFace,
 138                                const Handle(IntTools_Context)& theCtx);
 139
 140 static
 141   Standard_Real MaxDistance(const Handle(Geom_Curve)& theC,
 142                             const Standard_Real aT,
 143                             GeomAPI_ProjectPointOnSurf& theProjPS);
 144
 145 static
 146   Standard_Real FindMaxDistance(const Handle(Geom_Curve)& theC,
 147                                 const Standard_Real theFirst,
 148                                 const Standard_Real theLast,
 149                                 GeomAPI_ProjectPointOnSurf& theProjPS,
 150                                 const Standard_Real theEps);
 151
 152 static
 153   Standard_Real FindMaxDistance(const Handle(Geom_Curve)& theCurve,
 154                                 const Standard_Real theFirst,
 155                                 const Standard_Real theLast,
 156                                 const TopoDS_Face& theFace,
 157                                 const Handle(IntTools_Context)& theContext);
 158
 159 static
 160   void CorrectPlaneBoundaries(Standard_Real& aUmin,
 161                               Standard_Real& aUmax,
 162                               Standard_Real& aVmin,
 163                               Standard_Real& aVmax);
 164
 165 //=======================================================================
 166 //function :
 167 //purpose  :
 168 //=======================================================================
 169 IntTools_FaceFace::IntTools_FaceFace()
 170 {
 171   myIsDone=Standard_False;
 172   myTangentFaces=Standard_False;
 173   //
 174   myHS1 = new GeomAdaptor_HSurface ();
 175   myHS2 = new GeomAdaptor_HSurface ();
 176   myTolReached2d=0.;
 177   myTolReached3d=0.;
 178   myTolReal = 0.;
 179   myTolF1 = 0.;
 180   myTolF2 = 0.;
 181   myTol = 0.;
 182   myFuzzyValue = Precision::Confusion();
 183   SetParameters(Standard_True, Standard_True, Standard_True, 1.e-07);
 184 }
 185 //=======================================================================
 186 //function : SetContext
 187 //purpose  :
 188 //=======================================================================
 189 void IntTools_FaceFace::SetContext(const Handle(IntTools_Context)& aContext)
 190 {
 191   myContext=aContext;
 192 }
 193 //=======================================================================
 194 //function : Context
 195 //purpose  :
 196 //=======================================================================
 197 const Handle(IntTools_Context)& IntTools_FaceFace::Context()const
 198 {
 199   return myContext;
 200 }
 201 //=======================================================================
 202 //function : Face1
 203 //purpose  :
 204 //=======================================================================
 205 const TopoDS_Face& IntTools_FaceFace::Face1() const
 206 {
 207   return myFace1;
 208 }
 209 //=======================================================================
 210 //function : Face2
 211 //purpose  :
 212 //=======================================================================
 213 const TopoDS_Face& IntTools_FaceFace::Face2() const
 214 {
 215   return myFace2;
 216 }
 217 //=======================================================================
 218 //function : TangentFaces
 219 //purpose  :
 220 //=======================================================================
 221 Standard_Boolean IntTools_FaceFace::TangentFaces() const
 222 {
 223   return myTangentFaces;
 224 }
 225 //=======================================================================
 226 //function : Points
 227 //purpose  :
 228 //=======================================================================
 229 const IntTools_SequenceOfPntOn2Faces& IntTools_FaceFace::Points() const
 230 {
 231   return myPnts;
 232 }
 233 //=======================================================================
 234 //function : IsDone
 235 //purpose  :
 236 //=======================================================================
 237 Standard_Boolean IntTools_FaceFace::IsDone() const
 238 {
 239   return myIsDone;
 240 }
 241 //=======================================================================
 242 //function : TolReached3d
 243 //purpose  :
 244 //=======================================================================
 245 Standard_Real IntTools_FaceFace::TolReached3d() const
 246 {
 247   return myTolReached3d;
 248 }
 249 //=======================================================================
 250 //function : TolReal
 251 //purpose  :
 252 //=======================================================================
 253 Standard_Real IntTools_FaceFace::TolReal() const
 254 {
 255   return myTolReal;
 256 }
 257 //=======================================================================
 258 //function : Lines
 259 //purpose  : return lines of intersection
 260 //=======================================================================
 261 const IntTools_SequenceOfCurves& IntTools_FaceFace::Lines() const
 262 {
 263   StdFail_NotDone_Raise_if
 264     (!myIsDone,
 265      "IntTools_FaceFace::Lines() => myIntersector NOT DONE");
 266   return mySeqOfCurve;
 267 }
 268 //=======================================================================
 269 //function : TolReached2d
 270 //purpose  :
 271 //=======================================================================
 272 Standard_Real IntTools_FaceFace::TolReached2d() const
 273 {
 274   return myTolReached2d;
 275 }
 276 // =======================================================================
 277 // function: SetParameters
 278 //
 279 // =======================================================================
 280 void IntTools_FaceFace::SetParameters(const Standard_Boolean ToApproxC3d,
 281                                       const Standard_Boolean ToApproxC2dOnS1,
 282                                       const Standard_Boolean ToApproxC2dOnS2,
 283                                       const Standard_Real ApproximationTolerance)
 284 {
 285   myApprox = ToApproxC3d;
 286   myApprox1 = ToApproxC2dOnS1;
 287   myApprox2 = ToApproxC2dOnS2;
 288   myTolApprox = ApproximationTolerance;
 289 }
 290 //=======================================================================
 291 //function : SetFuzzyValue
 292 //purpose  :
 293 //=======================================================================
 294 void IntTools_FaceFace::SetFuzzyValue(const Standard_Real theFuzz)
 295 {
 296   myFuzzyValue = Max(theFuzz, Precision::Confusion());
 297 }
 298 //=======================================================================
 299 //function : FuzzyValue
 300 //purpose  :
 301 //=======================================================================
 302 Standard_Real IntTools_FaceFace::FuzzyValue() const
 303 {
 304   return myFuzzyValue;
 305 }
 306
 307 //=======================================================================
 308 //function : SetList
 309 //purpose  :
 310 //=======================================================================
 311 void IntTools_FaceFace::SetList(IntSurf_ListOfPntOn2S& aListOfPnts)
 312 {
 313   myListOfPnts = aListOfPnts;
 314 }
 315
 316
 317 static Standard_Boolean isTreatAnalityc(const BRepAdaptor_Surface& theBAS1,
 318                                         const BRepAdaptor_Surface& theBAS2,
 319                                         const Standard_Real theTol)
 320 {
 321   const Standard_Real Tolang = 1.e-8;
 322   Standard_Real aHigh = 0.0;
 323
 324   const GeomAbs_SurfaceType aType1=theBAS1.GetType();
 325   const GeomAbs_SurfaceType aType2=theBAS2.GetType();
 326
 327   gp_Pln aS1;
 328   gp_Cylinder aS2;
 329   if(aType1 == GeomAbs_Plane)
 330   {
 331     aS1=theBAS1.Plane();
 332   }
 333   else if(aType2 == GeomAbs_Plane)
 334   {
 335     aS1=theBAS2.Plane();
 336   }
 337   else
 338   {
 339     return Standard_True;
 340   }
 341
 342   if(aType1 == GeomAbs_Cylinder)
 343   {
 344     aS2=theBAS1.Cylinder();
 345     const Standard_Real VMin = theBAS1.FirstVParameter();
 346     const Standard_Real VMax = theBAS1.LastVParameter();
 347
 348     if( Precision::IsNegativeInfinite(VMin) ||
 349         Precision::IsPositiveInfinite(VMax))
 350           return Standard_True;
 351     else
 352       aHigh = VMax - VMin;
 353   }
 354   else if(aType2 == GeomAbs_Cylinder)
 355   {
 356     aS2=theBAS2.Cylinder();
 357
 358     const Standard_Real VMin = theBAS2.FirstVParameter();
 359     const Standard_Real VMax = theBAS2.LastVParameter();
 360
 361     if( Precision::IsNegativeInfinite(VMin) ||
 362         Precision::IsPositiveInfinite(VMax))
 363           return Standard_True;
 364     else
 365       aHigh = VMax - VMin;
 366   }
 367   else
 368   {
 369     return Standard_True;
 370   }
 371
 372   IntAna_QuadQuadGeo inter;
 373   inter.Perform(aS1,aS2,Tolang,theTol, aHigh);
 374   if(inter.TypeInter() == IntAna_Ellipse)
 375   {
 376     const gp_Elips anEl = inter.Ellipse(1);
 377     const Standard_Real aMajorR = anEl.MajorRadius();
 378     const Standard_Real aMinorR = anEl.MinorRadius();
 379
 380     return (aMajorR < 100000.0 * aMinorR);
 381   }
 382   else
 383   {
 384     return inter.IsDone();
 385   }
 386 }
 387 //=======================================================================
 388 //function : Perform
 389 //purpose  : intersect surfaces of the faces
 390 //=======================================================================
 391 void IntTools_FaceFace::Perform(const TopoDS_Face& aF1,
 392                                 const TopoDS_Face& aF2)
 393 {
 394   Standard_Boolean RestrictLine = Standard_False;
 395
 396   if (myContext.IsNull()) {
 397     myContext=new IntTools_Context;
 398   }
 399
 400   mySeqOfCurve.Clear();
 401   myTolReached2d=0.;
 402   myTolReached3d=0.;
 403   myTolReal = 0.;
 404   myIsDone = Standard_False;
 405   myNbrestr=0;//?
 406
 407   myFace1=aF1;
 408   myFace2=aF2;
 409
 410   const BRepAdaptor_Surface& aBAS1 = myContext->SurfaceAdaptor(myFace1);
 411   const BRepAdaptor_Surface& aBAS2 = myContext->SurfaceAdaptor(myFace2);
 412   GeomAbs_SurfaceType aType1=aBAS1.GetType();
 413   GeomAbs_SurfaceType aType2=aBAS2.GetType();
 414
 415   const Standard_Boolean bReverse=SortTypes(aType1, aType2);
 416   if (bReverse)
 417   {
 418     myFace1=aF2;
 419     myFace2=aF1;
 420     aType1=aBAS2.GetType();
 421     aType2=aBAS1.GetType();
 422
 423     if (myListOfPnts.Extent())
 424     {
 425       Standard_Real aU1,aV1,aU2,aV2;
 426       IntSurf_ListIteratorOfListOfPntOn2S aItP2S;
 427       //
 428       aItP2S.Initialize(myListOfPnts);
 429       for (; aItP2S.More(); aItP2S.Next())
 430       {
 431         IntSurf_PntOn2S& aP2S=aItP2S.Value();
 432         aP2S.Parameters(aU1,aV1,aU2,aV2);
 433         aP2S.SetValue(aU2,aV2,aU1,aV1);
 434       }
 435     }
 436     //
 437     Standard_Boolean anAproxTmp = myApprox1;
 438     myApprox1 = myApprox2;
 439     myApprox2 = anAproxTmp;
 440   }
 441
 442
 443   const Handle(Geom_Surface) S1=BRep_Tool::Surface(myFace1);
 444   const Handle(Geom_Surface) S2=BRep_Tool::Surface(myFace2);
 445
 446   Standard_Real aFuzz = myFuzzyValue / 2.;
 447   myTolF1 = BRep_Tool::Tolerance(myFace1) + aFuzz;
 448   myTolF2 = BRep_Tool::Tolerance(myFace2) + aFuzz;
 449   myTol = myTolF1 + myTolF2;
 450
 451   Standard_Real TolArc = myTol;
 452   Standard_Real TolTang = TolArc;
 453
 454   const Standard_Boolean isFace1Quad = (aType1 == GeomAbs_Cylinder ||
 455                                         aType1 == GeomAbs_Cone ||
 456                                         aType1 == GeomAbs_Torus);
 457
 458   const Standard_Boolean isFace2Quad = (aType2 == GeomAbs_Cylinder ||
 459                                         aType2 == GeomAbs_Cone ||
 460                                         aType2 == GeomAbs_Torus);
 461
 462   if(aType1==GeomAbs_Plane && aType2==GeomAbs_Plane)  {
 463     Standard_Real umin, umax, vmin, vmax;
 464     //
 465     myContext->UVBounds(myFace1, umin, umax, vmin, vmax);
 466     CorrectPlaneBoundaries(umin, umax, vmin, vmax);
 467     myHS1->ChangeSurface().Load(S1, umin, umax, vmin, vmax);
 468     //
 469     myContext->UVBounds(myFace2, umin, umax, vmin, vmax);
 470     CorrectPlaneBoundaries(umin, umax, vmin, vmax);
 471     myHS2->ChangeSurface().Load(S2, umin, umax, vmin, vmax);
 472     //
 473     Standard_Real TolAng = 1.e-8;
 474     //
 475     PerformPlanes(myHS1, myHS2,
 476                   myTolF1, myTolF2, TolAng, TolTang,
 477                   myApprox1, myApprox2,
 478                   mySeqOfCurve, myTangentFaces, myTolReached3d);
 479     //
 480     myIsDone = Standard_True;
 481
 482     if(!myTangentFaces) {
 483       const Standard_Integer NbLinPP = mySeqOfCurve.Length();
 484       if(NbLinPP) {
 485         Standard_Real aTolFMax;
 486         aTolFMax=Max(myTolF1, myTolF2);
 487         myTolReal = Precision::Confusion();
 488         if (aTolFMax > myTolReal) {
 489           myTolReal = aTolFMax;
 490         }
 491         if (aTolFMax > myTolReached3d) {
 492           myTolReached3d = aTolFMax;
 493         }
 494         //
 495         myTolReached2d = myTolReal;
 496
 497         if (bReverse) {
 498           Handle(Geom2d_Curve) aC2D1, aC2D2;
 499           const Standard_Integer aNbLin = mySeqOfCurve.Length();
 500           for (Standard_Integer i = 1; i <= aNbLin; ++i) {
 501             IntTools_Curve& aIC=mySeqOfCurve(i);
 502             aC2D1=aIC.FirstCurve2d();
 503             aC2D2=aIC.SecondCurve2d();
 504             aIC.SetFirstCurve2d(aC2D2);
 505             aIC.SetSecondCurve2d(aC2D1);
 506           }
 507         }
 508       }
 509     }
 510     return;
 511   }//if(aType1==GeomAbs_Plane && aType2==GeomAbs_Plane){
 512
 513   if ((aType1==GeomAbs_Plane) && isFace2Quad)
 514   {
 515     Standard_Real umin, umax, vmin, vmax;
 516     // F1
 517     myContext->UVBounds(myFace1, umin, umax, vmin, vmax);
 518     CorrectPlaneBoundaries(umin, umax, vmin, vmax);
 519     myHS1->ChangeSurface().Load(S1, umin, umax, vmin, vmax);
 520     // F2
 521     myContext->UVBounds(myFace2, umin, umax, vmin, vmax);
 522     CorrectSurfaceBoundaries(myFace2, myTol * 2., umin, umax, vmin, vmax);
 523     myHS2->ChangeSurface().Load(S2, umin, umax, vmin, vmax);
 524   }
 525   else if ((aType2==GeomAbs_Plane) && isFace1Quad)
 526   {
 527     Standard_Real umin, umax, vmin, vmax;
 528     //F1
 529     myContext->UVBounds(myFace1, umin, umax, vmin, vmax);
 530     CorrectSurfaceBoundaries(myFace1, myTol * 2., umin, umax, vmin, vmax);
 531     myHS1->ChangeSurface().Load(S1, umin, umax, vmin, vmax);
 532     // F2
 533     myContext->UVBounds(myFace2, umin, umax, vmin, vmax);
 534     CorrectPlaneBoundaries(umin, umax, vmin, vmax);
 535     myHS2->ChangeSurface().Load(S2, umin, umax, vmin, vmax);
 536   }
 537   else
 538   {
 539     Standard_Real umin, umax, vmin, vmax;
 540     myContext->UVBounds(myFace1, umin, umax, vmin, vmax);
 541     CorrectSurfaceBoundaries(myFace1, myTol * 2., umin, umax, vmin, vmax);
 542     myHS1->ChangeSurface().Load(S1, umin, umax, vmin, vmax);
 543     myContext->UVBounds(myFace2, umin, umax, vmin, vmax);
 544     CorrectSurfaceBoundaries(myFace2, myTol * 2., umin, umax, vmin, vmax);
 545     myHS2->ChangeSurface().Load(S2, umin, umax, vmin, vmax);
 546   }
 547
 548   const Handle(IntTools_TopolTool) dom1 = new IntTools_TopolTool(myHS1);
 549   const Handle(IntTools_TopolTool) dom2 = new IntTools_TopolTool(myHS2);
 550
 551   myLConstruct.Load(dom1, dom2, myHS1, myHS2);
 552
 553
 554   Tolerances(myHS1, myHS2, TolTang);
 555
 556   {
 557     const Standard_Real UVMaxStep = 0.001;
 558     const Standard_Real Deflection = 0.1;
 559     myIntersector.SetTolerances(TolArc, TolTang, UVMaxStep, Deflection);
 560   }
 561
 562   if((myHS1->IsUClosed() && !myHS1->IsUPeriodic()) ||
 563      (myHS1->IsVClosed() && !myHS1->IsVPeriodic()) ||
 564      (myHS2->IsUClosed() && !myHS2->IsUPeriodic()) ||
 565      (myHS2->IsVClosed() && !myHS2->IsVPeriodic()))
 566   {
 567     RestrictLine = Standard_True;
 568   }
 569   //
 570   if((aType1 != GeomAbs_BSplineSurface) &&
 571       (aType1 != GeomAbs_BezierSurface)  &&
 572      (aType1 != GeomAbs_OtherSurface)  &&
 573      (aType2 != GeomAbs_BSplineSurface) &&
 574       (aType2 != GeomAbs_BezierSurface)  &&
 575      (aType2 != GeomAbs_OtherSurface))
 576   {
 577     RestrictLine = Standard_True;
 578
 579     if ((aType1 == GeomAbs_Torus) ||
 580         (aType2 == GeomAbs_Torus))
 581     {
 582       myListOfPnts.Clear();
 583     }
 584   }
 585
 586   //
 587   if(!RestrictLine)
 588   {
 589     TopExp_Explorer aExp;
 590     for(Standard_Integer i = 0; (!RestrictLine) && (i < 2); i++)
 591     {
 592       const TopoDS_Face& aF=(!i) ? myFace1 : myFace2;
 593       aExp.Init(aF, TopAbs_EDGE);
 594       for(; aExp.More(); aExp.Next())
 595       {
 596         const TopoDS_Edge& aE=TopoDS::Edge(aExp.Current());
 597
 598         if(BRep_Tool::Degenerated(aE))
 599         {
 600           RestrictLine = Standard_True;
 601           break;
 602         }
 603       }
 604     }
 605   }
 606
 607 #ifdef INTTOOLS_FACEFACE_DEBUG
 608     if(!myListOfPnts.IsEmpty()) {
 609       char aBuff[10000];
 610       const IntSurf_PntOn2S& aPt = myListOfPnts.First();
 611       Standard_Real u1, v1, u2, v2;
 612       aPt.Parameters(u1, v1, u2, v2);
 613
 614       Sprintf(aBuff,"bopcurves <face1 face2> -2d");
 615       IntSurf_ListIteratorOfListOfPntOn2S IterLOP1(myListOfPnts);
 616       for(;IterLOP1.More(); IterLOP1.Next())
 617       {
 618         const IntSurf_PntOn2S& aPt = IterLOP1.Value();
 619         Standard_Real u1, v1, u2, v2;
 620         aPt.Parameters(u1, v1, u2, v2);
 621
 622         Sprintf(aBuff, "%s -p %+10.20f %+10.20f %+10.20f %+10.20f", aBuff, u1, v1, u2, v2);
 623       }
 624
 625       cout << aBuff << endl;
 626     }
 627 #endif
 628
 629   const Standard_Boolean isGeomInt = isTreatAnalityc(aBAS1, aBAS2, myTol);
 630   myIntersector.Perform(myHS1, dom1, myHS2, dom2, TolArc, TolTang,
 631                                   myListOfPnts, RestrictLine, isGeomInt);
 632
 633   myIsDone = myIntersector.IsDone();
 634
 635   if (myIsDone)
 636   {
 637     myTangentFaces=myIntersector.TangentFaces();
 638     if (myTangentFaces) {
 639       return;
 640     }
 641     //
 642     if(RestrictLine) {
 643       myListOfPnts.Clear(); // to use LineConstructor
 644     }
 645     //
 646     const Standard_Integer aNbLinIntersector = myIntersector.NbLines();
 647     for (Standard_Integer i=1; i <= aNbLinIntersector; ++i) {
 648       MakeCurve(i, dom1, dom2, TolArc);
 649     }
 650     //
 651     ComputeTolReached3d();
 652     //
 653     if (bReverse) {
 654       Handle(Geom2d_Curve) aC2D1, aC2D2;
 655       //
 656       const Standard_Integer aNbLinSeqOfCurve =mySeqOfCurve.Length();
 657       for (Standard_Integer i=1; i<=aNbLinSeqOfCurve; ++i)
 658       {
 659         IntTools_Curve& aIC=mySeqOfCurve(i);
 660         aC2D1=aIC.FirstCurve2d();
 661         aC2D2=aIC.SecondCurve2d();
 662         aIC.SetFirstCurve2d(aC2D2);
 663         aIC.SetSecondCurve2d(aC2D1);
 664       }
 665     }
 666
 667     // Points
 668     Standard_Boolean bValid2D1, bValid2D2;
 669     Standard_Real U1,V1,U2,V2;
 670     IntTools_PntOnFace aPntOnF1, aPntOnF2;
 671     IntTools_PntOn2Faces aPntOn2Faces;
 672     //
 673     const Standard_Integer aNbPnts = myIntersector.NbPnts();
 674     for (Standard_Integer i=1; i <= aNbPnts; ++i)
 675     {
 676       const IntSurf_PntOn2S& aISPnt=myIntersector.Point(i).PntOn2S();
 677       const gp_Pnt& aPnt=aISPnt.Value();
 678       aISPnt.Parameters(U1,V1,U2,V2);
 679       //
 680       // check the validity of the intersection point for the faces
 681       bValid2D1 = myContext->IsPointInOnFace(myFace1, gp_Pnt2d(U1, V1));
 682       if (!bValid2D1) {
 683         continue;
 684       }
 685       //
 686       bValid2D2 = myContext->IsPointInOnFace(myFace2, gp_Pnt2d(U2, V2));
 687       if (!bValid2D2) {
 688         continue;
 689       }
 690       //
 691       // add the intersection point
 692       aPntOnF1.Init(myFace1, aPnt, U1, V1);
 693       aPntOnF2.Init(myFace2, aPnt, U2, V2);
 694       //
 695       if (!bReverse)
 696       {
 697         aPntOn2Faces.SetP1(aPntOnF1);
 698         aPntOn2Faces.SetP2(aPntOnF2);
 699       }
 700       else
 701       {
 702         aPntOn2Faces.SetP2(aPntOnF1);
 703         aPntOn2Faces.SetP1(aPntOnF2);
 704       }
 705
 706       myPnts.Append(aPntOn2Faces);
 707     }
 708   }
 709 }
 710
 711 //=======================================================================
 712 //function : ComputeTolerance
 713 //purpose  :
 714 //=======================================================================
 715 Standard_Real IntTools_FaceFace::ComputeTolerance()
 716 {
 717   Standard_Integer i, j, aNbLin;
 718   Standard_Real aFirst, aLast, aD, aDMax, aT;
 719   Handle(Geom_Surface) aS1, aS2;
 720   //
 721   aDMax = 0;
 722   aNbLin = mySeqOfCurve.Length();
 723   //
 724   aS1 = myHS1->ChangeSurface().Surface();
 725   aS2 = myHS2->ChangeSurface().Surface();
 726   //
 727   for (i = 1; i <= aNbLin; ++i)
 728   {
 729     const IntTools_Curve& aIC = mySeqOfCurve(i);
 730     const Handle(Geom_Curve)& aC3D = aIC.Curve();
 731     if (aC3D.IsNull())
 732     {
 733       continue;
 734     }
 735     //
 736     aFirst = aC3D->FirstParameter();
 737     aLast  = aC3D->LastParameter();
 738     //
 739     const Handle(Geom2d_Curve)& aC2D1 = aIC.FirstCurve2d();
 740     const Handle(Geom2d_Curve)& aC2D2 = aIC.SecondCurve2d();
 741     //
 742     for (j = 0; j < 2; ++j)
 743     {
 744       const Handle(Geom2d_Curve)& aC2D = !j ? aC2D1 : aC2D2;
 745       const Handle(Geom_Surface)& aS = !j ? aS1 : aS2;
 746       //
 747       if (!aC2D.IsNull())
 748       {
 749         if (IntTools_Tools::ComputeTolerance
 750             (aC3D, aC2D, aS, aFirst, aLast, aD, aT))
 751         {
 752           if (aD > aDMax)
 753           {
 754             aDMax = aD;
 755           }
 756         }
 757       }
 758       else
 759       {
 760         const TopoDS_Face& aF = !j ? myFace1 : myFace2;
 761         aD = FindMaxDistance(aC3D, aFirst, aLast, aF, myContext);
 762         if (aD > aDMax)
 763         {
 764           aDMax = aD;
 765         }
 766       }
 767     }
 768   }
 769   //
 770   return aDMax;
 771 }
 772
 773 //=======================================================================
 774 //function :ComputeTolReached3d
 775 //purpose  :
 776 //=======================================================================
 777 void IntTools_FaceFace::ComputeTolReached3d()
 778 {
 779   Standard_Integer aNbLin;
 780   GeomAbs_SurfaceType aType1, aType2;
 781   //
 782   aNbLin=myIntersector.NbLines();
 783   if (!aNbLin) {
 784     return;
 785   }
 786   //
 787   aType1=myHS1->Surface().GetType();
 788   aType2=myHS2->Surface().GetType();
 789   //
 790   if (aType1==GeomAbs_Cylinder && aType2==GeomAbs_Cylinder)
 791   {
 792     if (aNbLin==2)
 793     {
 794       Handle(IntPatch_Line) aIL1, aIL2;
 795       IntPatch_IType aTL1, aTL2;
 796       //
 797       aIL1=myIntersector.Line(1);
 798       aIL2=myIntersector.Line(2);
 799       aTL1=aIL1->ArcType();
 800       aTL2=aIL2->ArcType();
 801       if (aTL1==IntPatch_Lin && aTL2==IntPatch_Lin) {
 802         Standard_Real aD, aDTresh, dTol;
 803         gp_Lin aL1, aL2;
 804         //
 805         dTol=1.e-8;
 806         aDTresh=1.5e-6;
 807         //
 808         aL1=Handle(IntPatch_GLine)::DownCast(aIL1)->Line();
 809         aL2=Handle(IntPatch_GLine)::DownCast(aIL2)->Line();
 810         aD=aL1.Distance(aL2);
 811         aD=0.5*aD;
 812         if (aD<aDTresh)
 813         {//In order to avoid creation too thin face
 814           myTolReached3d=aD+dTol;
 815         }
 816       }
 817     }
 818   }// if (aType1==GeomAbs_Cylinder && aType2==GeomAbs_Cylinder) {
 819   //
 820
 821   Standard_Real aDMax = ComputeTolerance();
 822   if (aDMax > myTolReached3d)
 823   {
 824     myTolReached3d = aDMax;
 825   }
 826   myTolReal = myTolReached3d;
 827 }
 828
 829 //=======================================================================
 830 //function : MakeCurve
 831 //purpose  :
 832 //=======================================================================
 833 void IntTools_FaceFace::MakeCurve(const Standard_Integer Index,
 834                                   const Handle(Adaptor3d_TopolTool)& dom1,
 835                                   const Handle(Adaptor3d_TopolTool)& dom2,
 836                                   const Standard_Real theToler)
 837 {
 838   Standard_Boolean bDone, rejectSurface, reApprox, bAvoidLineConstructor;
 839   Standard_Boolean ok, bPCurvesOk;
 840   Standard_Integer i, j, aNbParts;
 841   Standard_Real fprm, lprm;
 842   Standard_Real Tolpc;
 843   Handle(IntPatch_Line) L;
 844   IntPatch_IType typl;
 845   Handle(Geom_Curve) newc;
 846   //
 847   const Standard_Real TOLCHECK   =0.0000001;
 848   const Standard_Real TOLANGCHECK=0.1;
 849   //
 850   rejectSurface = Standard_False;
 851   reApprox = Standard_False;
 852   //
 853   bPCurvesOk = Standard_True;
 854
 855  reapprox:;
 856
 857   Tolpc = myTolApprox;
 858   bAvoidLineConstructor = Standard_False;
 859   L = myIntersector.Line(Index);
 860   typl = L->ArcType();
 861   //
 862   if(typl==IntPatch_Walking) {
 863     Handle(IntPatch_WLine) aWLine (Handle(IntPatch_WLine)::DownCast(L));
 864     if(aWLine.IsNull()) {
 865       return;
 866     }
 867     L = aWLine;
 868
 869     Standard_Integer nbp = aWLine->NbPnts();
 870     const IntSurf_PntOn2S& p1 = aWLine->Point(1);
 871     const IntSurf_PntOn2S& p2 = aWLine->Point(nbp);
 872
 873     const gp_Pnt& P1 = p1.Value();
 874     const gp_Pnt& P2 = p2.Value();
 875
 876     if(P1.SquareDistance(P2) < 1.e-14) {
 877       bAvoidLineConstructor = Standard_False;
 878     }
 879   }
 880
 881   typl=L->ArcType();
 882
 883   if(typl == IntPatch_Restriction)
 884     bAvoidLineConstructor = Standard_True;
 885
 886   //
 887   // Line Constructor
 888   if(!bAvoidLineConstructor) {
 889     myLConstruct.Perform(L);
 890     //
 891     bDone=myLConstruct.IsDone();
 892     if(!bDone)
 893     {
 894       return;
 895     }
 896
 897     if(typl != IntPatch_Restriction)
 898     {
 899       aNbParts=myLConstruct.NbParts();
 900       if (aNbParts <= 0)
 901       {
 902         return;
 903       }
 904     }
 905   }
 906   // Do the Curve
 907
 908
 909   switch (typl) {
 910   //########################################
 911   // Line, Parabola, Hyperbola
 912   //########################################
 913   case IntPatch_Lin:
 914   case IntPatch_Parabola:
 915   case IntPatch_Hyperbola: {
 916     if (typl == IntPatch_Lin) {
 917       newc =
 918         new Geom_Line (Handle(IntPatch_GLine)::DownCast(L)->Line());
 919     }
 920
 921     else if (typl == IntPatch_Parabola) {
 922       newc =
 923         new Geom_Parabola(Handle(IntPatch_GLine)::DownCast(L)->Parabola());
 924     }
 925
 926     else if (typl == IntPatch_Hyperbola) {
 927       newc =
 928         new Geom_Hyperbola (Handle(IntPatch_GLine)::DownCast(L)->Hyperbola());
 929     }
 930     //
 931     // myTolReached3d
 932     if (typl == IntPatch_Lin) {
 933       TolR3d (myTolF1, myTolF2, myTolReached3d);
 934     }
 935     //
 936     aNbParts=myLConstruct.NbParts();
 937     for (i=1; i<=aNbParts; i++) {
 938       Standard_Boolean bFNIt, bLPIt;
 939       //
 940       myLConstruct.Part(i, fprm, lprm);
 941         //
 942       bFNIt=Precision::IsNegativeInfinite(fprm);
 943       bLPIt=Precision::IsPositiveInfinite(lprm);
 944       //
 945       if (!bFNIt && !bLPIt) {
 946         //
 947         IntTools_Curve aCurve;
 948         //
 949         Handle(Geom_TrimmedCurve) aCT3D=new Geom_TrimmedCurve(newc, fprm, lprm);
 950         aCurve.SetCurve(aCT3D);
 951         if (typl == IntPatch_Parabola) {
 952           myTolReached3d=IntTools_Tools::CurveTolerance(aCT3D, myTol);
 953         }
 954         //
 955         aCurve.SetCurve(new Geom_TrimmedCurve(newc, fprm, lprm));
 956         if(myApprox1) {
 957           Handle (Geom2d_Curve) C2d;
 958           GeomInt_IntSS::BuildPCurves(fprm, lprm, Tolpc,
 959                 myHS1->ChangeSurface().Surface(), newc, C2d);
 960           if(Tolpc>myTolReached2d || myTolReached2d==0.) {
 961             myTolReached2d=Tolpc;
 962           }
 963             //
 964             aCurve.SetFirstCurve2d(new Geom2d_TrimmedCurve(C2d,fprm,lprm));
 965           }
 966           else {
 967             Handle(Geom2d_BSplineCurve) H1;
 968             //
 969             aCurve.SetFirstCurve2d(H1);
 970           }
 971         //
 972         if(myApprox2) {
 973           Handle (Geom2d_Curve) C2d;
 974           GeomInt_IntSS::BuildPCurves(fprm, lprm, Tolpc,
 975                     myHS2->ChangeSurface().Surface(), newc, C2d);
 976           if(Tolpc>myTolReached2d || myTolReached2d==0.) {
 977             myTolReached2d=Tolpc;
 978           }
 979           //
 980           aCurve.SetSecondCurve2d(new Geom2d_TrimmedCurve(C2d,fprm,lprm));
 981           }
 982         else {
 983           Handle(Geom2d_BSplineCurve) H1;
 984           //
 985           aCurve.SetSecondCurve2d(H1);
 986         }
 987         mySeqOfCurve.Append(aCurve);
 988       } //if (!bFNIt && !bLPIt) {
 989       else {
 990         //  on regarde si on garde
 991         //
 992         Standard_Real aTestPrm, dT=100.;
 993         //
 994         aTestPrm=0.;
 995         if (bFNIt && !bLPIt) {
 996           aTestPrm=lprm-dT;
 997         }
 998         else if (!bFNIt && bLPIt) {
 999           aTestPrm=fprm+dT;
1000         }
1001         else {
1002           // i.e, if (bFNIt && bLPIt)
1003           aTestPrm=IntTools_Tools::IntermediatePoint(-dT, dT);
1004         }
1005         //
1006         gp_Pnt ptref(newc->Value(aTestPrm));
1007         //
1008         GeomAbs_SurfaceType typS1 = myHS1->GetType();
1009         GeomAbs_SurfaceType typS2 = myHS2->GetType();
1010         if( typS1 == GeomAbs_SurfaceOfExtrusion ||
1011             typS1 == GeomAbs_OffsetSurface ||
1012             typS1 == GeomAbs_SurfaceOfRevolution ||
1013             typS2 == GeomAbs_SurfaceOfExtrusion ||
1014             typS2 == GeomAbs_OffsetSurface ||
1015             typS2 == GeomAbs_SurfaceOfRevolution) {
1016           Handle(Geom2d_BSplineCurve) H1;
1017           mySeqOfCurve.Append(IntTools_Curve(newc, H1, H1));
1018           continue;
1019         }
1020
1021         Standard_Real u1, v1, u2, v2, Tol;
1022
1023         Tol = Precision::Confusion();
1024         Parameters(myHS1, myHS2, ptref,  u1, v1, u2, v2);
1025         ok = (dom1->Classify(gp_Pnt2d(u1, v1), Tol) != TopAbs_OUT);
1026         if(ok) {
1027           ok = (dom2->Classify(gp_Pnt2d(u2,v2),Tol) != TopAbs_OUT);
1028         }
1029         if (ok) {
1030           Handle(Geom2d_BSplineCurve) H1;
1031           mySeqOfCurve.Append(IntTools_Curve(newc, H1, H1));
1032         }
1033       }
1034     }// for (i=1; i<=aNbParts; i++) {
1035   }// case IntPatch_Lin:  case IntPatch_Parabola:  case IntPatch_Hyperbola:
1036     break;
1037
1038   //########################################
1039   // Circle and Ellipse
1040   //########################################
1041   case IntPatch_Circle:
1042   case IntPatch_Ellipse: {
1043
1044     if (typl == IntPatch_Circle) {
1045       newc = new Geom_Circle
1046         (Handle(IntPatch_GLine)::DownCast(L)->Circle());
1047     }
1048     else { //IntPatch_Ellipse
1049       newc = new Geom_Ellipse
1050         (Handle(IntPatch_GLine)::DownCast(L)->Ellipse());
1051     }
1052     //
1053     // myTolReached3d
1054     TolR3d (myTolF1, myTolF2, myTolReached3d);
1055     //
1056     aNbParts=myLConstruct.NbParts();
1057     //
1058     Standard_Real aPeriod, aNul;
1059     TColStd_SequenceOfReal aSeqFprm,  aSeqLprm;
1060
1061     aNul=0.;
1062     aPeriod=M_PI+M_PI;
1063
1064     for (i=1; i<=aNbParts; i++) {
1065       myLConstruct.Part(i, fprm, lprm);
1066
1067       if (fprm < aNul && lprm > aNul) {
1068         // interval that goes through 0. is divided on two intervals;
1069         while (fprm<aNul || fprm>aPeriod) fprm=fprm+aPeriod;
1070         while (lprm<aNul || lprm>aPeriod) lprm=lprm+aPeriod;
1071         //
1072         if((aPeriod - fprm) > Tolpc) {
1073           aSeqFprm.Append(fprm);
1074           aSeqLprm.Append(aPeriod);
1075         }
1076         else {
1077           gp_Pnt P1 = newc->Value(fprm);
1078           gp_Pnt P2 = newc->Value(aPeriod);
1079           Standard_Real aTolDist = myTol;
1080           aTolDist = (myTolReached3d > aTolDist) ? myTolReached3d : aTolDist;
1081
1082           if(P1.Distance(P2) > aTolDist) {
1083             Standard_Real anewpar = fprm;
1084
1085             if(ParameterOutOfBoundary(fprm, newc, myFace1, myFace2,
1086                                       lprm, Standard_False, myTol, anewpar, myContext)) {
1087               fprm = anewpar;
1088             }
1089             aSeqFprm.Append(fprm);
1090             aSeqLprm.Append(aPeriod);
1091           }
1092         }
1093
1094         //
1095         if((lprm - aNul) > Tolpc) {
1096           aSeqFprm.Append(aNul);
1097           aSeqLprm.Append(lprm);
1098         }
1099         else {
1100           gp_Pnt P1 = newc->Value(aNul);
1101           gp_Pnt P2 = newc->Value(lprm);
1102           Standard_Real aTolDist = myTol;
1103           aTolDist = (myTolReached3d > aTolDist) ? myTolReached3d : aTolDist;
1104
1105           if(P1.Distance(P2) > aTolDist) {
1106             Standard_Real anewpar = lprm;
1107
1108             if(ParameterOutOfBoundary(lprm, newc, myFace1, myFace2,
1109                                       fprm, Standard_True, myTol, anewpar, myContext)) {
1110               lprm = anewpar;
1111             }
1112             aSeqFprm.Append(aNul);
1113             aSeqLprm.Append(lprm);
1114           }
1115         }
1116       }
1117       else {
1118         // usual interval
1119         aSeqFprm.Append(fprm);
1120         aSeqLprm.Append(lprm);
1121       }
1122     }
1123
1124     //
1125     aNbParts=aSeqFprm.Length();
1126     for (i=1; i<=aNbParts; i++) {
1127       fprm=aSeqFprm(i);
1128       lprm=aSeqLprm(i);
1129       //
1130       Standard_Real aRealEpsilon=RealEpsilon();
1131       if (Abs(fprm) > aRealEpsilon || Abs(lprm-2.*M_PI) > aRealEpsilon) {
1132         //==============================================
1133         ////
1134         IntTools_Curve aCurve;
1135         Handle(Geom_TrimmedCurve) aTC3D=new Geom_TrimmedCurve(newc,fprm,lprm);
1136         aCurve.SetCurve(aTC3D);
1137         fprm=aTC3D->FirstParameter();
1138         lprm=aTC3D->LastParameter ();
1139         ////
1140         if (typl == IntPatch_Circle || typl == IntPatch_Ellipse) {////
1141           if(myApprox1) {
1142             Handle (Geom2d_Curve) C2d;
1143             GeomInt_IntSS::BuildPCurves(fprm, lprm, Tolpc,
1144                         myHS1->ChangeSurface().Surface(), newc, C2d);
1145             if(Tolpc>myTolReached2d || myTolReached2d==0) {
1146               myTolReached2d=Tolpc;
1147             }
1148             //
1149             aCurve.SetFirstCurve2d(C2d);
1150           }
1151           else { ////
1152             Handle(Geom2d_BSplineCurve) H1;
1153             aCurve.SetFirstCurve2d(H1);
1154           }
1155
1156
1157           if(myApprox2) {
1158             Handle (Geom2d_Curve) C2d;
1159             GeomInt_IntSS::BuildPCurves(fprm,lprm,Tolpc,
1160                     myHS2->ChangeSurface().Surface(),newc,C2d);
1161             if(Tolpc>myTolReached2d || myTolReached2d==0) {
1162               myTolReached2d=Tolpc;
1163             }
1164             //
1165             aCurve.SetSecondCurve2d(C2d);
1166           }
1167           else {
1168             Handle(Geom2d_BSplineCurve) H1;
1169             aCurve.SetSecondCurve2d(H1);
1170           }
1171         }
1172
1173         else {
1174           Handle(Geom2d_BSplineCurve) H1;
1175           aCurve.SetFirstCurve2d(H1);
1176           aCurve.SetSecondCurve2d(H1);
1177         }
1178         mySeqOfCurve.Append(aCurve);
1179           //==============================================
1180       } //if (Abs(fprm) > RealEpsilon() || Abs(lprm-2.*M_PI) > RealEpsilon())
1181
1182       else {
1183         //  on regarde si on garde
1184         //
1185         if (aNbParts==1) {
1186 //           if (Abs(fprm) < RealEpsilon() &&  Abs(lprm-2.*M_PI) < RealEpsilon()) {
1187           if (Abs(fprm) <= aRealEpsilon && Abs(lprm-2.*M_PI) <= aRealEpsilon) {
1188             IntTools_Curve aCurve;
1189             Handle(Geom_TrimmedCurve) aTC3D=new Geom_TrimmedCurve(newc,fprm,lprm);
1190             aCurve.SetCurve(aTC3D);
1191             fprm=aTC3D->FirstParameter();
1192             lprm=aTC3D->LastParameter ();
1193
1194             if(myApprox1) {
1195               Handle (Geom2d_Curve) C2d;
1196               GeomInt_IntSS::BuildPCurves(fprm,lprm,Tolpc,
1197                     myHS1->ChangeSurface().Surface(),newc,C2d);
1198               if(Tolpc>myTolReached2d || myTolReached2d==0) {
1199                 myTolReached2d=Tolpc;
1200               }
1201               //
1202               aCurve.SetFirstCurve2d(C2d);
1203             }
1204             else { ////
1205               Handle(Geom2d_BSplineCurve) H1;
1206               aCurve.SetFirstCurve2d(H1);
1207             }
1208
1209             if(myApprox2) {
1210               Handle (Geom2d_Curve) C2d;
1211               GeomInt_IntSS::BuildPCurves(fprm,lprm,Tolpc,
1212                     myHS2->ChangeSurface().Surface(),newc,C2d);
1213               if(Tolpc>myTolReached2d || myTolReached2d==0) {
1214                 myTolReached2d=Tolpc;
1215               }
1216               //
1217               aCurve.SetSecondCurve2d(C2d);
1218             }
1219             else {
1220               Handle(Geom2d_BSplineCurve) H1;
1221               aCurve.SetSecondCurve2d(H1);
1222             }
1223             mySeqOfCurve.Append(aCurve);
1224             break;
1225           }
1226         }
1227         //
1228         Standard_Real aTwoPIdiv17, u1, v1, u2, v2, Tol;
1229
1230         aTwoPIdiv17=2.*M_PI/17.;
1231
1232         for (j=0; j<=17; j++) {
1233           gp_Pnt ptref (newc->Value (j*aTwoPIdiv17));
1234           Tol = Precision::Confusion();
1235
1236           Parameters(myHS1, myHS2, ptref, u1, v1, u2, v2);
1237           ok = (dom1->Classify(gp_Pnt2d(u1,v1),Tol) != TopAbs_OUT);
1238           if(ok) {
1239             ok = (dom2->Classify(gp_Pnt2d(u2,v2),Tol) != TopAbs_OUT);
1240           }
1241           if (ok) {
1242             IntTools_Curve aCurve;
1243             aCurve.SetCurve(newc);
1244             //==============================================
1245             if (typl == IntPatch_Circle || typl == IntPatch_Ellipse) {
1246
1247               if(myApprox1) {
1248                 Handle (Geom2d_Curve) C2d;
1249                 GeomInt_IntSS::BuildPCurves(fprm, lprm, Tolpc,
1250                         myHS1->ChangeSurface().Surface(), newc, C2d);
1251                 if(Tolpc>myTolReached2d || myTolReached2d==0) {
1252                   myTolReached2d=Tolpc;
1253                 }
1254                 //
1255                 aCurve.SetFirstCurve2d(C2d);
1256               }
1257               else {
1258                 Handle(Geom2d_BSplineCurve) H1;
1259                 aCurve.SetFirstCurve2d(H1);
1260               }
1261
1262               if(myApprox2) {
1263                 Handle (Geom2d_Curve) C2d;
1264                 GeomInt_IntSS::BuildPCurves(fprm, lprm, Tolpc,
1265                         myHS2->ChangeSurface().Surface(), newc, C2d);
1266                 if(Tolpc>myTolReached2d || myTolReached2d==0) {
1267                   myTolReached2d=Tolpc;
1268                 }
1269                 //
1270                 aCurve.SetSecondCurve2d(C2d);
1271               }
1272
1273               else {
1274                 Handle(Geom2d_BSplineCurve) H1;
1275                 aCurve.SetSecondCurve2d(H1);
1276               }
1277             }//  end of if (typl == IntPatch_Circle || typl == IntPatch_Ellipse)
1278
1279             else {
1280               Handle(Geom2d_BSplineCurve) H1;
1281               //
1282               aCurve.SetFirstCurve2d(H1);
1283               aCurve.SetSecondCurve2d(H1);
1284             }
1285             //==============================================
1286             //
1287             mySeqOfCurve.Append(aCurve);
1288             break;
1289
1290             }//  end of if (ok) {
1291           }//  end of for (Standard_Integer j=0; j<=17; j++)
1292         }//  end of else { on regarde si on garde
1293       }// for (i=1; i<=myLConstruct.NbParts(); i++)
1294     }// IntPatch_Circle: IntPatch_Ellipse:
1295     break;
1296
1297   case IntPatch_Analytic:
1298     //This case was processed earlier (in IntPatch_Intersection)
1299     break;
1300
1301   case IntPatch_Walking:{
1302     Handle(IntPatch_WLine) WL =
1303       Handle(IntPatch_WLine)::DownCast(L);
1304
1305 #ifdef INTTOOLS_FACEFACE_DEBUG
1306     WL->Dump(0);
1307 #endif
1308
1309     //
1310     Standard_Integer ifprm, ilprm;
1311     //
1312     if (!myApprox) {
1313       aNbParts = 1;
1314       if(!bAvoidLineConstructor){
1315         aNbParts=myLConstruct.NbParts();
1316       }
1317       for (i=1; i<=aNbParts; ++i) {
1318         Handle(Geom2d_BSplineCurve) H1, H2;
1319         Handle(Geom_Curve) aBSp;
1320         //
1321         if(bAvoidLineConstructor) {
1322           ifprm = 1;
1323           ilprm = WL->NbPnts();
1324         }
1325         else {
1326           myLConstruct.Part(i, fprm, lprm);
1327           ifprm=(Standard_Integer)fprm;
1328           ilprm=(Standard_Integer)lprm;
1329         }
1330         //
1331         if(myApprox1) {
1332           H1 = GeomInt_IntSS::MakeBSpline2d(WL, ifprm, ilprm, Standard_True);
1333         }
1334         //
1335         if(myApprox2) {
1336           H2 = GeomInt_IntSS::MakeBSpline2d(WL, ifprm, ilprm, Standard_False);
1337         }
1338         //
1339         aBSp=GeomInt_IntSS::MakeBSpline(WL, ifprm, ilprm);
1340         IntTools_Curve aIC(aBSp, H1, H2);
1341         mySeqOfCurve.Append(aIC);
1342       }// for (i=1; i<=aNbParts; ++i) {
1343     }// if (!myApprox) {
1344     //
1345     else { // X
1346       Standard_Boolean bIsDecomposited;
1347       Standard_Integer nbiter, aNbSeqOfL;
1348       Standard_Real tol2d, aTolApproxImp;
1349       IntPatch_SequenceOfLine aSeqOfL;
1350       GeomInt_WLApprox theapp3d;
1351       Approx_ParametrizationType aParType = Approx_ChordLength;
1352       //
1353       Standard_Boolean anApprox1 = myApprox1;
1354       Standard_Boolean anApprox2 = myApprox2;
1355       //
1356       aTolApproxImp=1.e-5;
1357       tol2d = myTolApprox;
1358
1359       GeomAbs_SurfaceType typs1, typs2;
1360       typs1 = myHS1->Surface().GetType();
1361       typs2 = myHS2->Surface().GetType();
1362       Standard_Boolean anWithPC = Standard_True;
1363
1364       if(typs1 == GeomAbs_Cylinder && typs2 == GeomAbs_Sphere) {
1365         anWithPC =
1366           ApproxWithPCurves(myHS1->Surface().Cylinder(), myHS2->Surface().Sphere());
1367       }
1368       else if (typs1 == GeomAbs_Sphere && typs2 == GeomAbs_Cylinder) {
1369         anWithPC =
1370           ApproxWithPCurves(myHS2->Surface().Cylinder(), myHS1->Surface().Sphere());
1371       }
1372       //
1373       if(!anWithPC) {
1374         myTolApprox = aTolApproxImp;//1.e-5;
1375         anApprox1 = Standard_False;
1376         anApprox2 = Standard_False;
1377         //
1378         tol2d = myTolApprox;
1379       }
1380
1381       if(myHS1 == myHS2) {
1382         theapp3d.SetParameters(myTolApprox, tol2d, 4, 8, 0, 30, Standard_False, aParType);
1383         rejectSurface = Standard_True;
1384       }
1385       else {
1386         if(reApprox && !rejectSurface)
1387           theapp3d.SetParameters(myTolApprox, tol2d, 4, 8, 0, 30, Standard_False, aParType);
1388         else {
1389           Standard_Integer iDegMax, iDegMin, iNbIter;
1390           //
1391           ApproxParameters(myHS1, myHS2, iDegMin, iDegMax, iNbIter);
1392           theapp3d.SetParameters(myTolApprox, tol2d, iDegMin, iDegMax,
1393                                                   iNbIter, 30, Standard_True, aParType);
1394         }
1395       }
1396       //
1397       Standard_Real aReachedTol = Precision::Confusion();
1398       bIsDecomposited = IntTools_WLineTool::
1399         DecompositionOfWLine(WL,
1400                              myHS1,
1401                              myHS2,
1402                              myFace1,
1403                              myFace2,
1404                              myLConstruct,
1405                              bAvoidLineConstructor,
1406                              myTol,
1407                              aSeqOfL,
1408                              aReachedTol,
1409                              myContext);
1410       if ( bIsDecomposited && ( myTolReached3d < aReachedTol ) ) {
1411         myTolReached3d = aReachedTol;
1412       }
1413       //
1414       aNbSeqOfL=aSeqOfL.Length();
1415       //
1416       if (bIsDecomposited) {
1417         nbiter=aNbSeqOfL;
1418       }
1419       else {
1420         nbiter=1;
1421         aNbParts=1;
1422         if (!bAvoidLineConstructor) {
1423           aNbParts=myLConstruct.NbParts();
1424           nbiter=aNbParts;
1425         }
1426       }
1427       //
1428       for(i = 1; i <= nbiter; ++i) {
1429         if(bIsDecomposited) {
1430           WL = Handle(IntPatch_WLine)::DownCast(aSeqOfL.Value(i));
1431           ifprm = 1;
1432           ilprm = WL->NbPnts();
1433         }
1434         else {
1435           if(bAvoidLineConstructor) {
1436             ifprm = 1;
1437             ilprm = WL->NbPnts();
1438           }
1439           else {
1440             myLConstruct.Part(i, fprm, lprm);
1441             ifprm = (Standard_Integer)fprm;
1442             ilprm = (Standard_Integer)lprm;
1443           }
1444         }
1445         //-- lbr :
1446         //-- Si une des surfaces est un plan , on approxime en 2d
1447         //-- sur cette surface et on remonte les points 2d en 3d.
1448         if(typs1 == GeomAbs_Plane) {
1449           theapp3d.Perform(myHS1, myHS2, WL, Standard_False,Standard_True, myApprox2,ifprm,ilprm);
1450         }
1451         else if(typs2 == GeomAbs_Plane) {
1452           theapp3d.Perform(myHS1,myHS2,WL,Standard_False,myApprox1,Standard_True,ifprm,ilprm);
1453         }
1454         else {
1455           //
1456           if (myHS1 != myHS2){
1457             if ((typs1==GeomAbs_BezierSurface || typs1==GeomAbs_BSplineSurface) &&
1458                 (typs2==GeomAbs_BezierSurface || typs2==GeomAbs_BSplineSurface)) {
1459
1460               theapp3d.SetParameters(myTolApprox, tol2d, 4, 8, 0, 30,
1461                                                                 Standard_True, aParType);
1462
1463               Standard_Boolean bUseSurfaces;
1464               bUseSurfaces = IntTools_WLineTool::NotUseSurfacesForApprox(myFace1, myFace2, WL, ifprm,  ilprm);
1465               if (bUseSurfaces) {
1466                 // ######
1467                 rejectSurface = Standard_True;
1468                 // ######
1469                 theapp3d.SetParameters(myTolApprox, tol2d, 4, 8, 0, 30,
1470                                                                 Standard_False, aParType);
1471               }
1472             }
1473           }
1474           //
1475           theapp3d.Perform(myHS1,myHS2,WL,Standard_True,anApprox1,anApprox2,ifprm,ilprm);
1476         }
1477           //
1478         if (!theapp3d.IsDone()) {
1479           Handle(Geom2d_BSplineCurve) H1;
1480           Handle(Geom2d_BSplineCurve) H2;
1481           //
1482           Handle(Geom_Curve) aBSp=GeomInt_IntSS::MakeBSpline(WL,ifprm, ilprm);
1483           //
1484           if(myApprox1) {
1485             H1 = GeomInt_IntSS::MakeBSpline2d(WL, ifprm, ilprm, Standard_True);
1486           }
1487           //
1488           if(myApprox2) {
1489             H2 = GeomInt_IntSS::MakeBSpline2d(WL, ifprm, ilprm, Standard_False);
1490           }
1491           //
1492           IntTools_Curve aIC(aBSp, H1, H2);
1493           mySeqOfCurve.Append(aIC);
1494         }
1495
1496         else {
1497           if(myApprox1 || myApprox2 || (typs1==GeomAbs_Plane || typs2==GeomAbs_Plane)) {
1498             if( theapp3d.TolReached2d()>myTolReached2d || myTolReached2d==0.) {
1499               myTolReached2d = theapp3d.TolReached2d();
1500             }
1501           }
1502           if(typs1==GeomAbs_Plane || typs2==GeomAbs_Plane) {
1503             myTolReached3d = myTolReached2d;
1504             //
1505             if (typs1==GeomAbs_Torus || typs2==GeomAbs_Torus) {
1506               if (myTolReached3d<1.e-6) {
1507                 myTolReached3d = theapp3d.TolReached3d();
1508                 myTolReached3d=1.e-6;
1509               }
1510             }
1511           }
1512           else  if( theapp3d.TolReached3d()>myTolReached3d || myTolReached3d==0.) {
1513             myTolReached3d = theapp3d.TolReached3d();
1514           }
1515
1516           Standard_Integer aNbMultiCurves, nbpoles;
1517           aNbMultiCurves=theapp3d.NbMultiCurves();
1518           for (j=1; j<=aNbMultiCurves; j++) {
1519             if(typs1 == GeomAbs_Plane) {
1520               const AppParCurves_MultiBSpCurve& mbspc = theapp3d.Value(j);
1521               nbpoles = mbspc.NbPoles();
1522
1523               TColgp_Array1OfPnt2d tpoles2d(1,nbpoles);
1524               TColgp_Array1OfPnt   tpoles(1,nbpoles);
1525
1526               mbspc.Curve(1,tpoles2d);
1527               const gp_Pln&  Pln = myHS1->Surface().Plane();
1528               //
1529               Standard_Integer ik;
1530               for(ik = 1; ik<= nbpoles; ik++) {
1531                 tpoles.SetValue(ik,
1532                                 ElSLib::Value(tpoles2d.Value(ik).X(),
1533                                               tpoles2d.Value(ik).Y(),
1534                                               Pln));
1535               }
1536               //
1537               Handle(Geom_BSplineCurve) BS =
1538                 new Geom_BSplineCurve(tpoles,
1539                                       mbspc.Knots(),
1540                                       mbspc.Multiplicities(),
1541                                       mbspc.Degree());
1542               GeomLib_CheckBSplineCurve Check(BS,TOLCHECK,TOLANGCHECK);
1543               Check.FixTangent(Standard_True, Standard_True);
1544               //
1545               IntTools_Curve aCurve;
1546               aCurve.SetCurve(BS);
1547
1548               if(myApprox1) {
1549                 Handle(Geom2d_BSplineCurve) BS1 =
1550                   new Geom2d_BSplineCurve(tpoles2d,
1551                                           mbspc.Knots(),
1552                                           mbspc.Multiplicities(),
1553                                           mbspc.Degree());
1554                 GeomLib_Check2dBSplineCurve Check1(BS1,TOLCHECK,TOLANGCHECK);
1555                 Check1.FixTangent(Standard_True,Standard_True);
1556                 //
1557                 // ############################################
1558                 if(!rejectSurface && !reApprox) {
1559                   Standard_Boolean isValid = IsCurveValid(BS1);
1560                   if(!isValid) {
1561                     reApprox = Standard_True;
1562                     goto reapprox;
1563                   }
1564                 }
1565                 // ############################################
1566                 aCurve.SetFirstCurve2d(BS1);
1567               }
1568               else {
1569                 Handle(Geom2d_BSplineCurve) H1;
1570                 aCurve.SetFirstCurve2d(H1);
1571               }
1572
1573               if(myApprox2) {
1574                 mbspc.Curve(2, tpoles2d);
1575
1576                 Handle(Geom2d_BSplineCurve) BS2 = new Geom2d_BSplineCurve(tpoles2d,
1577                                                                           mbspc.Knots(),
1578                                                                           mbspc.Multiplicities(),
1579                                                                           mbspc.Degree());
1580                 GeomLib_Check2dBSplineCurve newCheck(BS2,TOLCHECK,TOLANGCHECK);
1581                 newCheck.FixTangent(Standard_True,Standard_True);
1582
1583                 // ###########################################
1584                 if(!rejectSurface && !reApprox) {
1585                   Standard_Boolean isValid = IsCurveValid(BS2);
1586                   if(!isValid) {
1587                     reApprox = Standard_True;
1588                     goto reapprox;
1589                   }
1590                 }
1591                 // ###########################################
1592                 //
1593                 aCurve.SetSecondCurve2d(BS2);
1594               }
1595               else {
1596                 Handle(Geom2d_BSplineCurve) H2;
1597                 //
1598                   aCurve.SetSecondCurve2d(H2);
1599               }
1600               //
1601               mySeqOfCurve.Append(aCurve);
1602
1603             }//if(typs1 == GeomAbs_Plane) {
1604
1605             else if(typs2 == GeomAbs_Plane)
1606             {
1607               const AppParCurves_MultiBSpCurve& mbspc = theapp3d.Value(j);
1608               nbpoles = mbspc.NbPoles();
1609
1610               TColgp_Array1OfPnt2d tpoles2d(1,nbpoles);
1611               TColgp_Array1OfPnt   tpoles(1,nbpoles);
1612               mbspc.Curve((myApprox1==Standard_True)? 2 : 1,tpoles2d);
1613               const gp_Pln&  Pln = myHS2->Surface().Plane();
1614               //
1615               Standard_Integer ik;
1616               for(ik = 1; ik<= nbpoles; ik++) {
1617                 tpoles.SetValue(ik,
1618                                 ElSLib::Value(tpoles2d.Value(ik).X(),
1619                                               tpoles2d.Value(ik).Y(),
1620                                               Pln));
1621
1622               }
1623               //
1624               Handle(Geom_BSplineCurve) BS=new Geom_BSplineCurve(tpoles,
1625                                                                  mbspc.Knots(),
1626                                                                  mbspc.Multiplicities(),
1627                                                                  mbspc.Degree());
1628               GeomLib_CheckBSplineCurve Check(BS,TOLCHECK,TOLANGCHECK);
1629               Check.FixTangent(Standard_True,Standard_True);
1630               //
1631               IntTools_Curve aCurve;
1632               aCurve.SetCurve(BS);
1633
1634               if(myApprox2) {
1635                 Handle(Geom2d_BSplineCurve) BS1=new Geom2d_BSplineCurve(tpoles2d,
1636                                                                         mbspc.Knots(),
1637                                                                         mbspc.Multiplicities(),
1638                                                                         mbspc.Degree());
1639                 GeomLib_Check2dBSplineCurve Check1(BS1,TOLCHECK,TOLANGCHECK);
1640                 Check1.FixTangent(Standard_True,Standard_True);
1641                 //
1642                 // ###########################################
1643                 if(!rejectSurface && !reApprox) {
1644                   Standard_Boolean isValid = IsCurveValid(BS1);
1645                   if(!isValid) {
1646                     reApprox = Standard_True;
1647                     goto reapprox;
1648                   }
1649                 }
1650                 // ###########################################
1651                 bPCurvesOk = CheckPCurve(BS1, myFace2, myContext);
1652                 aCurve.SetSecondCurve2d(BS1);
1653               }
1654               else {
1655                 Handle(Geom2d_BSplineCurve) H2;
1656                 aCurve.SetSecondCurve2d(H2);
1657               }
1658
1659               if(myApprox1) {
1660                 mbspc.Curve(1,tpoles2d);
1661                 Handle(Geom2d_BSplineCurve) BS2=new Geom2d_BSplineCurve(tpoles2d,
1662                                                                         mbspc.Knots(),
1663                                                                         mbspc.Multiplicities(),
1664                                                                         mbspc.Degree());
1665                 GeomLib_Check2dBSplineCurve Check2(BS2,TOLCHECK,TOLANGCHECK);
1666                 Check2.FixTangent(Standard_True,Standard_True);
1667                 //
1668                 // ###########################################
1669                 if(!rejectSurface && !reApprox) {
1670                   Standard_Boolean isValid = IsCurveValid(BS2);
1671                   if(!isValid) {
1672                     reApprox = Standard_True;
1673                     goto reapprox;
1674                   }
1675                 }
1676                 // ###########################################
1677                 bPCurvesOk = bPCurvesOk && CheckPCurve(BS2, myFace1, myContext);
1678                 aCurve.SetFirstCurve2d(BS2);
1679               }
1680               else {
1681                 Handle(Geom2d_BSplineCurve) H1;
1682                 //
1683                 aCurve.SetFirstCurve2d(H1);
1684               }
1685               //
1686               //if points of the pcurves are out of the faces bounds
1687               //create 3d and 2d curves without approximation
1688               if (!bPCurvesOk) {
1689                 Handle(Geom2d_BSplineCurve) H1, H2;
1690                 bPCurvesOk = Standard_True;
1691                 //
1692                 Handle(Geom_Curve) aBSp=GeomInt_IntSS::MakeBSpline(WL,ifprm, ilprm);
1693
1694                 if(myApprox1) {
1695                   H1 = GeomInt_IntSS::MakeBSpline2d(WL, ifprm, ilprm, Standard_True);
1696                   bPCurvesOk = CheckPCurve(H1, myFace1, myContext);
1697                 }
1698
1699                 if(myApprox2) {
1700                   H2 = GeomInt_IntSS::MakeBSpline2d(WL, ifprm, ilprm, Standard_False);
1701                   bPCurvesOk = bPCurvesOk && CheckPCurve(H2, myFace2, myContext);
1702                 }
1703                 //
1704                 //if pcurves created without approximation are out of the
1705                 //faces bounds, use approximated 3d and 2d curves
1706                 if (bPCurvesOk) {
1707                   IntTools_Curve aIC(aBSp, H1, H2);
1708                   mySeqOfCurve.Append(aIC);
1709                 } else {
1710                   mySeqOfCurve.Append(aCurve);
1711                 }
1712               } else {
1713                 mySeqOfCurve.Append(aCurve);
1714               }
1715
1716             }// else if(typs2 == GeomAbs_Plane)
1717             //
1718             else { //typs2 != GeomAbs_Plane && typs1 != GeomAbs_Plane
1719               Standard_Boolean bIsValid1, bIsValid2;
1720               Handle(Geom_BSplineCurve) BS;
1721               Handle(Geom2d_BSplineCurve) aH2D;
1722               IntTools_Curve aCurve;
1723               //
1724               bIsValid1=Standard_True;
1725               bIsValid2=Standard_True;
1726               //
1727               const AppParCurves_MultiBSpCurve& mbspc = theapp3d.Value(j);
1728               nbpoles = mbspc.NbPoles();
1729               TColgp_Array1OfPnt tpoles(1,nbpoles);
1730               mbspc.Curve(1,tpoles);
1731               BS=new Geom_BSplineCurve(tpoles,
1732                                                                  mbspc.Knots(),
1733                                                                  mbspc.Multiplicities(),
1734                                                                  mbspc.Degree());
1735               GeomLib_CheckBSplineCurve Check(BS,TOLCHECK,TOLANGCHECK);
1736               Check.FixTangent(Standard_True,Standard_True);
1737               //
1738               aCurve.SetCurve(BS);
1739               aCurve.SetFirstCurve2d(aH2D);
1740               aCurve.SetSecondCurve2d(aH2D);
1741               //
1742               if(myApprox1) {
1743                 if(anApprox1) {
1744                   Handle(Geom2d_BSplineCurve) BS1;
1745                   TColgp_Array1OfPnt2d tpoles2d(1,nbpoles);
1746                   mbspc.Curve(2,tpoles2d);
1747                   //
1748                   BS1=new Geom2d_BSplineCurve(tpoles2d,
1749                                                                         mbspc.Knots(),
1750                                                                         mbspc.Multiplicities(),
1751                                                                         mbspc.Degree());
1752                   GeomLib_Check2dBSplineCurve newCheck(BS1,TOLCHECK,TOLANGCHECK);
1753                   newCheck.FixTangent(Standard_True,Standard_True);
1754                   //
1755                   if (!reApprox) {
1756                     bIsValid1=CheckPCurve(BS1, myFace1, myContext);
1757                   }
1758                   //
1759                   aCurve.SetFirstCurve2d(BS1);
1760                 }
1761                 else {
1762                   Handle(Geom2d_BSplineCurve) BS1;
1763                   fprm = BS->FirstParameter();
1764                   lprm = BS->LastParameter();
1765
1766                   Handle(Geom2d_Curve) C2d;
1767                   Standard_Real aTol = myTolApprox;
1768                   GeomInt_IntSS::BuildPCurves(fprm, lprm, aTol,
1769                             myHS1->ChangeSurface().Surface(), BS, C2d);
1770                   BS1 = Handle(Geom2d_BSplineCurve)::DownCast(C2d);
1771                   aCurve.SetFirstCurve2d(BS1);
1772                 }
1773               } // if(myApprox1) {
1774                 //
1775               if(myApprox2) {
1776                 if(anApprox2) {
1777                   Handle(Geom2d_BSplineCurve) BS2;
1778                   TColgp_Array1OfPnt2d tpoles2d(1,nbpoles);
1779                   mbspc.Curve((myApprox1==Standard_True)? 3 : 2,tpoles2d);
1780                   BS2=new Geom2d_BSplineCurve(tpoles2d,
1781                                                                         mbspc.Knots(),
1782                                                                         mbspc.Multiplicities(),
1783                                                                         mbspc.Degree());
1784                   GeomLib_Check2dBSplineCurve newCheck(BS2,TOLCHECK,TOLANGCHECK);
1785                   newCheck.FixTangent(Standard_True,Standard_True);
1786                 //
1787                   if (!reApprox) {
1788                     bIsValid2=CheckPCurve(BS2, myFace2, myContext);
1789                   }
1790                   aCurve.SetSecondCurve2d(BS2);
1791                 }
1792                 else {
1793                   Handle(Geom2d_BSplineCurve) BS2;
1794                   fprm = BS->FirstParameter();
1795                   lprm = BS->LastParameter();
1796
1797                   Handle(Geom2d_Curve) C2d;
1798                   Standard_Real aTol = myTolApprox;
1799                   GeomInt_IntSS::BuildPCurves(fprm, lprm, aTol,
1800                             myHS2->ChangeSurface().Surface(), BS, C2d);
1801                   BS2 = Handle(Geom2d_BSplineCurve)::DownCast(C2d);
1802                   aCurve.SetSecondCurve2d(BS2);
1803                 }
1804               } //if(myApprox2) {
1805               if (!bIsValid1 || !bIsValid2) {
1806                 myTolApprox=aTolApproxImp;//1.e-5;
1807                 tol2d = myTolApprox;
1808                 reApprox = Standard_True;
1809                 goto reapprox;
1810               }
1811                 //
1812               mySeqOfCurve.Append(aCurve);
1813             }
1814           }
1815         }
1816       }
1817     }// else { // X
1818   }// case IntPatch_Walking:{
1819     break;
1820
1821   case IntPatch_Restriction:
1822     {
1823       Handle(IntPatch_RLine) RL =
1824         Handle(IntPatch_RLine)::DownCast(L);
1825
1826 #ifdef INTTOOLS_FACEFACE_DEBUG
1827     RL->Dump(0);
1828 #endif
1829
1830       Handle(Geom_Curve) aC3d;
1831       Handle(Geom2d_Curve) aC2d1, aC2d2;
1832       Standard_Real aTolReached;
1833       GeomInt_IntSS::TreatRLine(RL, myHS1, myHS2, aC3d,
1834                                   aC2d1, aC2d2, aTolReached);
1835
1836       if(aC3d.IsNull())
1837         break;
1838
1839       Bnd_Box2d aBox1, aBox2;
1840
1841       const Standard_Real aU1f = myHS1->FirstUParameter(),
1842                           aV1f = myHS1->FirstVParameter(),
1843                           aU1l = myHS1->LastUParameter(),
1844                           aV1l = myHS1->LastVParameter();
1845       const Standard_Real aU2f = myHS2->FirstUParameter(),
1846                           aV2f = myHS2->FirstVParameter(),
1847                           aU2l = myHS2->LastUParameter(),
1848                           aV2l = myHS2->LastVParameter();
1849
1850       aBox1.Add(gp_Pnt2d(aU1f, aV1f));
1851       aBox1.Add(gp_Pnt2d(aU1l, aV1l));
1852       aBox2.Add(gp_Pnt2d(aU2f, aV2f));
1853       aBox2.Add(gp_Pnt2d(aU2l, aV2l));
1854
1855       GeomInt_VectorOfReal anArrayOfParameters;
1856
1857       //We consider here that the intersection line is same-parameter-line
1858       anArrayOfParameters.Append(aC3d->FirstParameter());
1859       anArrayOfParameters.Append(aC3d->LastParameter());
1860
1861       GeomInt_IntSS::
1862         TrimILineOnSurfBoundaries(aC2d1, aC2d2, aBox1, aBox2, anArrayOfParameters);
1863
1864       //Intersect with true boundaries. After that, enlarge bounding-boxes in order to
1865       //correct definition, if point on curve is inscribed in the box.
1866       aBox1.Enlarge(theToler);
1867       aBox2.Enlarge(theToler);
1868
1869       const Standard_Integer aNbIntersSolutionsm1 = anArrayOfParameters.Length() - 1;
1870
1871       //Trim RLine found.
1872       for(Standard_Integer anInd = 0; anInd < aNbIntersSolutionsm1; anInd++)
1873       {
1874         Standard_Real &aParF = anArrayOfParameters(anInd),
1875                       &aParL = anArrayOfParameters(anInd+1);
1876
1877         if((aParL - aParF) <= Precision::PConfusion())
1878         {
1879           //In order to more precise extending to the boundaries of source curves.
1880           if(anInd < aNbIntersSolutionsm1-1)
1881             aParL = aParF;
1882
1883           continue;
1884         }
1885
1886         const Standard_Real aPar = 0.5*(aParF + aParL);
1887         gp_Pnt2d aPt;
1888
1889         Handle(Geom2d_Curve) aCurv2d1, aCurv2d2;
1890         if(!aC2d1.IsNull())
1891         {
1892           aC2d1->D0(aPar, aPt);
1893
1894           if(aBox1.IsOut(aPt))
1895             continue;
1896
1897           if(myApprox1)
1898             aCurv2d1 = new Geom2d_TrimmedCurve(aC2d1, aParF, aParL);
1899         }
1900
1901         if(!aC2d2.IsNull())
1902         {
1903           aC2d2->D0(aPar, aPt);
1904
1905           if(aBox2.IsOut(aPt))
1906             continue;
1907
1908           if(myApprox2)
1909             aCurv2d2 = new Geom2d_TrimmedCurve(aC2d2, aParF, aParL);
1910         }
1911
1912         Handle(Geom_Curve) aCurv3d = new Geom_TrimmedCurve(aC3d, aParF, aParL);
1913
1914         IntTools_Curve aIC(aCurv3d, aCurv2d1, aCurv2d2);
1915         mySeqOfCurve.Append(aIC);
1916       }
1917     }
1918     break;
1919   default:
1920     break;
1921
1922   }
1923 }
1924
1925 //=======================================================================
1926 //function : Parameters
1927 //purpose  :
1928 //=======================================================================
1929  void Parameters(const Handle(GeomAdaptor_HSurface)& HS1,
1930                  const Handle(GeomAdaptor_HSurface)& HS2,
1931                  const gp_Pnt& Ptref,
1932                  Standard_Real& U1,
1933                  Standard_Real& V1,
1934                  Standard_Real& U2,
1935                  Standard_Real& V2)
1936 {
1937
1938   IntSurf_Quadric quad1,quad2;
1939   GeomAbs_SurfaceType typs = HS1->Surface().GetType();
1940
1941   switch (typs) {
1942   case GeomAbs_Plane:
1943     quad1.SetValue(HS1->Surface().Plane());
1944     break;
1945   case GeomAbs_Cylinder:
1946     quad1.SetValue(HS1->Surface().Cylinder());
1947     break;
1948   case GeomAbs_Cone:
1949     quad1.SetValue(HS1->Surface().Cone());
1950     break;
1951   case GeomAbs_Sphere:
1952     quad1.SetValue(HS1->Surface().Sphere());
1953     break;
1954   case GeomAbs_Torus:
1955     quad1.SetValue(HS1->Surface().Torus());
1956     break;
1957   default:
1958     Standard_ConstructionError::Raise("GeomInt_IntSS::MakeCurve");
1959   }
1960
1961   typs = HS2->Surface().GetType();
1962   switch (typs) {
1963   case GeomAbs_Plane:
1964     quad2.SetValue(HS2->Surface().Plane());
1965     break;
1966   case GeomAbs_Cylinder:
1967     quad2.SetValue(HS2->Surface().Cylinder());
1968     break;
1969   case GeomAbs_Cone:
1970     quad2.SetValue(HS2->Surface().Cone());
1971     break;
1972   case GeomAbs_Sphere:
1973     quad2.SetValue(HS2->Surface().Sphere());
1974     break;
1975   case GeomAbs_Torus:
1976     quad2.SetValue(HS2->Surface().Torus());
1977     break;
1978   default:
1979     Standard_ConstructionError::Raise("GeomInt_IntSS::MakeCurve");
1980   }
1981
1982   quad1.Parameters(Ptref,U1,V1);
1983   quad2.Parameters(Ptref,U2,V2);
1984 }
1985
1986 //=======================================================================
1987 //function : MakeBSpline
1988 //purpose  :
1989 //=======================================================================
1990 Handle(Geom_Curve) MakeBSpline  (const Handle(IntPatch_WLine)& WL,
1991                                  const Standard_Integer ideb,
1992                                  const Standard_Integer ifin)
1993 {
1994   Standard_Integer i,nbpnt = ifin-ideb+1;
1995   TColgp_Array1OfPnt poles(1,nbpnt);
1996   TColStd_Array1OfReal knots(1,nbpnt);
1997   TColStd_Array1OfInteger mults(1,nbpnt);
1998   Standard_Integer ipidebm1;
1999   for(i=1,ipidebm1=i+ideb-1; i<=nbpnt;ipidebm1++, i++) {
2000     poles(i) = WL->Point(ipidebm1).Value();
2001     mults(i) = 1;
2002     knots(i) = i-1;
2003   }
2004   mults(1) = mults(nbpnt) = 2;
2005   return
2006     new Geom_BSplineCurve(poles,knots,mults,1);
2007 }
2008
2009 //=======================================================================
2010 //function : PrepareLines3D
2011 //purpose  :
2012 //=======================================================================
2013   void IntTools_FaceFace::PrepareLines3D(const Standard_Boolean bToSplit)
2014 {
2015   Standard_Integer i, aNbCurves;
2016   GeomAbs_SurfaceType aType1, aType2;
2017   IntTools_SequenceOfCurves aNewCvs;
2018   //
2019   // 1. Treatment closed  curves
2020   aNbCurves=mySeqOfCurve.Length();
2021   for (i=1; i<=aNbCurves; ++i) {
2022     const IntTools_Curve& aIC=mySeqOfCurve(i);
2023     //
2024     if (bToSplit) {
2025       Standard_Integer j, aNbC;
2026       IntTools_SequenceOfCurves aSeqCvs;
2027       //
2028       aNbC=IntTools_Tools::SplitCurve(aIC, aSeqCvs);
2029       if (aNbC) {
2030         for (j=1; j<=aNbC; ++j) {
2031           const IntTools_Curve& aICNew=aSeqCvs(j);
2032           aNewCvs.Append(aICNew);
2033         }
2034       }
2035       else {
2036         aNewCvs.Append(aIC);
2037       }
2038     }
2039     else {
2040       aNewCvs.Append(aIC);
2041     }
2042   }
2043   //
2044   // 2. Plane\Cone intersection when we had 4 curves
2045   aType1=myHS1->GetType();
2046   aType2=myHS2->GetType();
2047   aNbCurves=aNewCvs.Length();
2048   //
2049   if ((aType1==GeomAbs_Plane && aType2==GeomAbs_Cone) ||
2050       (aType2==GeomAbs_Plane && aType1==GeomAbs_Cone)) {
2051     if (aNbCurves==4) {
2052       GeomAbs_CurveType aCType1;
2053       //
2054       aCType1=aNewCvs(1).Type();
2055       if (aCType1==GeomAbs_Line) {
2056         IntTools_SequenceOfCurves aSeqIn, aSeqOut;
2057         //
2058         for (i=1; i<=aNbCurves; ++i) {
2059           const IntTools_Curve& aIC=aNewCvs(i);
2060           aSeqIn.Append(aIC);
2061         }
2062         //
2063         IntTools_Tools::RejectLines(aSeqIn, aSeqOut);
2064         //
2065         aNewCvs.Clear();
2066         aNbCurves=aSeqOut.Length();
2067         for (i=1; i<=aNbCurves; ++i) {
2068           const IntTools_Curve& aIC=aSeqOut(i);
2069           aNewCvs.Append(aIC);
2070         }
2071       }
2072     }
2073   }// if ((aType1==GeomAbs_Plane && aType2==GeomAbs_Cone)...
2074   //
2075   // 3. Fill  mySeqOfCurve
2076   mySeqOfCurve.Clear();
2077   aNbCurves=aNewCvs.Length();
2078   for (i=1; i<=aNbCurves; ++i) {
2079     const IntTools_Curve& aIC=aNewCvs(i);
2080     mySeqOfCurve.Append(aIC);
2081   }
2082 }
2083 //=======================================================================
2084 //function : CorrectSurfaceBoundaries
2085 //purpose  :
2086 //=======================================================================
2087  void CorrectSurfaceBoundaries(const TopoDS_Face&  theFace,
2088                               const Standard_Real theTolerance,
2089                               Standard_Real&      theumin,
2090                               Standard_Real&      theumax,
2091                               Standard_Real&      thevmin,
2092                               Standard_Real&      thevmax)
2093 {
2094   Standard_Boolean enlarge, isuperiodic, isvperiodic;
2095   Standard_Real uinf, usup, vinf, vsup, delta;
2096   GeomAbs_SurfaceType aType;
2097   Handle(Geom_Surface) aSurface;
2098   //
2099   aSurface = BRep_Tool::Surface(theFace);
2100   aSurface->Bounds(uinf, usup, vinf, vsup);
2101   delta = theTolerance;
2102   enlarge = Standard_False;
2103   //
2104   GeomAdaptor_Surface anAdaptorSurface(aSurface);
2105   //
2106   if(aSurface->IsKind(STANDARD_TYPE(Geom_RectangularTrimmedSurface))) {
2107     Handle(Geom_Surface) aBasisSurface =
2108       (Handle(Geom_RectangularTrimmedSurface)::DownCast(aSurface))->BasisSurface();
2109
2110     if(aBasisSurface->IsKind(STANDARD_TYPE(Geom_RectangularTrimmedSurface)) ||
2111        aBasisSurface->IsKind(STANDARD_TYPE(Geom_OffsetSurface))) {
2112       return;
2113     }
2114   }
2115   //
2116   if(aSurface->IsKind(STANDARD_TYPE(Geom_OffsetSurface))) {
2117     Handle(Geom_Surface) aBasisSurface =
2118       (Handle(Geom_OffsetSurface)::DownCast(aSurface))->BasisSurface();
2119
2120     if(aBasisSurface->IsKind(STANDARD_TYPE(Geom_RectangularTrimmedSurface)) ||
2121        aBasisSurface->IsKind(STANDARD_TYPE(Geom_OffsetSurface))) {
2122       return;
2123     }
2124   }
2125   //
2126   isuperiodic = anAdaptorSurface.IsUPeriodic();
2127   isvperiodic = anAdaptorSurface.IsVPeriodic();
2128   //
2129   aType=anAdaptorSurface.GetType();
2130   if((aType==GeomAbs_BezierSurface) ||
2131      (aType==GeomAbs_BSplineSurface) ||
2132      (aType==GeomAbs_SurfaceOfExtrusion) ||
2133      (aType==GeomAbs_SurfaceOfRevolution) ||
2134      (aType==GeomAbs_Cylinder)) {
2135     enlarge=Standard_True;
2136   }
2137   //
2138   if(!isuperiodic && enlarge) {
2139
2140     if(!Precision::IsInfinite(theumin) &&
2141         ((theumin - uinf) > delta))
2142       theumin -= delta;
2143     else {
2144       theumin = uinf;
2145     }
2146
2147     if(!Precision::IsInfinite(theumax) &&
2148         ((usup - theumax) > delta))
2149       theumax += delta;
2150     else
2151       theumax = usup;
2152   }
2153   //
2154   if(!isvperiodic && enlarge) {
2155     if(!Precision::IsInfinite(thevmin) &&
2156         ((thevmin - vinf) > delta)) {
2157       thevmin -= delta;
2158     }
2159     else {
2160       thevmin = vinf;
2161     }
2162     if(!Precision::IsInfinite(thevmax) &&
2163         ((vsup - thevmax) > delta)) {
2164       thevmax += delta;
2165     }
2166     else {
2167       thevmax = vsup;
2168     }
2169   }
2170   //
2171   if(isuperiodic || isvperiodic) {
2172     Standard_Boolean correct = Standard_False;
2173     Standard_Boolean correctU = Standard_False;
2174     Standard_Boolean correctV = Standard_False;
2175     Bnd_Box2d aBox;
2176     TopExp_Explorer anExp;
2177
2178     for(anExp.Init(theFace, TopAbs_EDGE); anExp.More(); anExp.Next()) {
2179       if(BRep_Tool::IsClosed(TopoDS::Edge(anExp.Current()), theFace)) {
2180         correct = Standard_True;
2181         Standard_Real f, l;
2182         TopoDS_Edge anEdge = TopoDS::Edge(anExp.Current());
2183
2184         for(Standard_Integer i = 0; i < 2; i++) {
2185           if(i==0) {
2186             anEdge.Orientation(TopAbs_FORWARD);
2187           }
2188           else {
2189             anEdge.Orientation(TopAbs_REVERSED);
2190           }
2191           Handle(Geom2d_Curve) aCurve = BRep_Tool::CurveOnSurface(anEdge, theFace, f, l);
2192
2193           if(aCurve.IsNull()) {
2194             correct = Standard_False;
2195             break;
2196           }
2197           Handle(Geom2d_Line) aLine = Handle(Geom2d_Line)::DownCast(aCurve);
2198
2199           if(aLine.IsNull()) {
2200             correct = Standard_False;
2201             break;
2202           }
2203           gp_Dir2d anUDir(1., 0.);
2204           gp_Dir2d aVDir(0., 1.);
2205           Standard_Real anAngularTolerance = Precision::Angular();
2206
2207           correctU = correctU || aLine->Position().Direction().IsParallel(aVDir, anAngularTolerance);
2208           correctV = correctV || aLine->Position().Direction().IsParallel(anUDir, anAngularTolerance);
2209
2210           gp_Pnt2d pp1 = aCurve->Value(f);
2211           aBox.Add(pp1);
2212           gp_Pnt2d pp2 = aCurve->Value(l);
2213           aBox.Add(pp2);
2214         }
2215         if(!correct)
2216           break;
2217       }
2218     }
2219
2220     if(correct) {
2221       Standard_Real umin, vmin, umax, vmax;
2222       aBox.Get(umin, vmin, umax, vmax);
2223
2224       if(isuperiodic && correctU) {
2225         if(theumin < umin)
2226           theumin = umin;
2227         if(theumax > umax) {
2228           theumax = umax;
2229         }
2230       }
2231       if(isvperiodic && correctV) {
2232         if(thevmin < vmin)
2233           thevmin = vmin;
2234         if(thevmax > vmax)
2235           thevmax = vmax;
2236       }
2237     }
2238   }
2239 }
2240
2241 //=======================================================================
2242 //function : TolR3d
2243 //purpose  :
2244 //=======================================================================
2245 void TolR3d(const Standard_Real aTolF1,
2246             const Standard_Real aTolF2,
2247             Standard_Real& myTolReached3d)
2248 {
2249   Standard_Real aTolFMax, aTolTresh;
2250
2251   aTolTresh=2.999999e-3;
2252   aTolFMax=Max(aTolF1, aTolF2);
2253
2254   if (aTolFMax>aTolTresh) {
2255     myTolReached3d=aTolFMax;
2256   }
2257 }
2258
2259 // ------------------------------------------------------------------------------------------------
2260 // static function: ParameterOutOfBoundary
2261 // purpose:         Computes a new parameter for given curve. The corresponding 2d points
2262 //                  does not lay on any boundary of given faces
2263 // ------------------------------------------------------------------------------------------------
2264 Standard_Boolean ParameterOutOfBoundary(const Standard_Real       theParameter,
2265                                         const Handle(Geom_Curve)& theCurve,
2266                                         const TopoDS_Face&        theFace1,
2267                                         const TopoDS_Face&        theFace2,
2268                                         const Standard_Real       theOtherParameter,
2269                                         const Standard_Boolean    bIncreasePar,
2270                                         const Standard_Real       theTol,
2271                                         Standard_Real&            theNewParameter,
2272                                         const Handle(IntTools_Context)& aContext)
2273 {
2274   Standard_Boolean bIsComputed = Standard_False;
2275   theNewParameter = theParameter;
2276
2277   Standard_Real acurpar = theParameter;
2278   TopAbs_State aState = TopAbs_ON;
2279   Standard_Integer iter = 0;
2280   Standard_Real asumtol = theTol;
2281   Standard_Real adelta = asumtol * 0.1;
2282   adelta = (adelta < Precision::Confusion()) ? Precision::Confusion() : adelta;
2283   Handle(Geom_Surface) aSurf1 = BRep_Tool::Surface(theFace1);
2284   Handle(Geom_Surface) aSurf2 = BRep_Tool::Surface(theFace2);
2285
2286   Standard_Real u1, u2, v1, v2;
2287
2288   GeomAPI_ProjectPointOnSurf aPrj1;
2289   aSurf1->Bounds(u1, u2, v1, v2);
2290   aPrj1.Init(aSurf1, u1, u2, v1, v2);
2291
2292   GeomAPI_ProjectPointOnSurf aPrj2;
2293   aSurf2->Bounds(u1, u2, v1, v2);
2294   aPrj2.Init(aSurf2, u1, u2, v1, v2);
2295
2296   while(aState == TopAbs_ON) {
2297     if(bIncreasePar)
2298       acurpar += adelta;
2299     else
2300       acurpar -= adelta;
2301     gp_Pnt aPCurrent = theCurve->Value(acurpar);
2302     aPrj1.Perform(aPCurrent);
2303     Standard_Real U=0., V=0.;
2304
2305     if(aPrj1.IsDone()) {
2306       aPrj1.LowerDistanceParameters(U, V);
2307       aState = aContext->StatePointFace(theFace1, gp_Pnt2d(U, V));
2308     }
2309
2310     if(aState != TopAbs_ON) {
2311       aPrj2.Perform(aPCurrent);
2312
2313       if(aPrj2.IsDone()) {
2314         aPrj2.LowerDistanceParameters(U, V);
2315         aState = aContext->StatePointFace(theFace2, gp_Pnt2d(U, V));
2316       }
2317     }
2318
2319     if(iter > 11) {
2320       break;
2321     }
2322     iter++;
2323   }
2324
2325   if(iter <= 11) {
2326     theNewParameter = acurpar;
2327     bIsComputed = Standard_True;
2328
2329     if(bIncreasePar) {
2330       if(acurpar >= theOtherParameter)
2331         theNewParameter = theOtherParameter;
2332     }
2333     else {
2334       if(acurpar <= theOtherParameter)
2335         theNewParameter = theOtherParameter;
2336     }
2337   }
2338   return bIsComputed;
2339 }
2340
2341 //=======================================================================
2342 //function : IsCurveValid
2343 //purpose  :
2344 //=======================================================================
2345 Standard_Boolean IsCurveValid (const Handle(Geom2d_Curve)& thePCurve)
2346 {
2347   if(thePCurve.IsNull())
2348     return Standard_False;
2349
2350   Standard_Real tolint = 1.e-10;
2351   Geom2dAdaptor_Curve PCA;
2352   IntRes2d_Domain PCD;
2353   Geom2dInt_GInter PCI;
2354
2355   Standard_Real pf = 0., pl = 0.;
2356   gp_Pnt2d pntf, pntl;
2357
2358   if(!thePCurve->IsClosed() && !thePCurve->IsPeriodic()) {
2359     pf = thePCurve->FirstParameter();
2360     pl = thePCurve->LastParameter();
2361     pntf = thePCurve->Value(pf);
2362     pntl = thePCurve->Value(pl);
2363     PCA.Load(thePCurve);
2364     if(!PCA.IsPeriodic()) {
2365       if(PCA.FirstParameter() > pf) pf = PCA.FirstParameter();
2366       if(PCA.LastParameter()  < pl) pl = PCA.LastParameter();
2367     }
2368     PCD.SetValues(pntf,pf,tolint,pntl,pl,tolint);
2369     PCI.Perform(PCA,PCD,tolint,tolint);
2370     if(PCI.IsDone())
2371       if(PCI.NbPoints() > 0) {
2372         return Standard_False;
2373       }
2374   }
2375
2376   return Standard_True;
2377 }
2378
2379 //=======================================================================
2380 //static function : ApproxWithPCurves
2381 //purpose  : for bug 20964 only
2382 //=======================================================================
2383 Standard_Boolean ApproxWithPCurves(const gp_Cylinder& theCyl,
2384                                    const gp_Sphere& theSph)
2385 {
2386   Standard_Boolean bRes = Standard_True;
2387   Standard_Real R1 = theCyl.Radius(), R2 = theSph.Radius();
2388   //
2389   {
2390     Standard_Real aD2, aRc2, aEps;
2391     gp_Pnt aApexSph;
2392     //
2393     aEps=1.E-7;
2394     aRc2=R1*R1;
2395     //
2396     const gp_Ax3& aAx3Sph=theSph.Position();
2397     const gp_Pnt& aLocSph=aAx3Sph.Location();
2398     const gp_Dir& aDirSph=aAx3Sph.Direction();
2399     //
2400     const gp_Ax1& aAx1Cyl=theCyl.Axis();
2401     gp_Lin aLinCyl(aAx1Cyl);
2402     //
2403     aApexSph.SetXYZ(aLocSph.XYZ()+R2*aDirSph.XYZ());
2404     aD2=aLinCyl.SquareDistance(aApexSph);
2405     if (fabs(aD2-aRc2)<aEps) {
2406       return !bRes;
2407     }
2408     //
2409     aApexSph.SetXYZ(aLocSph.XYZ()-R2*aDirSph.XYZ());
2410     aD2=aLinCyl.SquareDistance(aApexSph);
2411     if (fabs(aD2-aRc2)<aEps) {
2412       return !bRes;
2413     }
2414   }
2415   //
2416
2417   if(R1 < 2.*R2) {
2418     return bRes;
2419   }
2420   gp_Lin anCylAx(theCyl.Axis());
2421
2422   Standard_Real aDist = anCylAx.Distance(theSph.Location());
2423   Standard_Real aDRel = Abs(aDist - R1)/R2;
2424
2425   if(aDRel > .2) return bRes;
2426
2427   Standard_Real par = ElCLib::Parameter(anCylAx, theSph.Location());
2428   gp_Pnt aP = ElCLib::Value(par, anCylAx);
2429   gp_Vec aV(aP, theSph.Location());
2430
2431   Standard_Real dd = aV.Dot(theSph.Position().XDirection());
2432
2433   if(aDist < R1 && dd > 0.) return Standard_False;
2434   if(aDist > R1 && dd < 0.) return Standard_False;
2435
2436
2437   return bRes;
2438 }
2439 //=======================================================================
2440 //function : PerformPlanes
2441 //purpose  :
2442 //=======================================================================
2443 void  PerformPlanes(const Handle(GeomAdaptor_HSurface)& theS1,
2444                     const Handle(GeomAdaptor_HSurface)& theS2,
2445                     const Standard_Real TolF1,
2446                     const Standard_Real TolF2,
2447                     const Standard_Real TolAng,
2448                     const Standard_Real TolTang,
2449                     const Standard_Boolean theApprox1,
2450                     const Standard_Boolean theApprox2,
2451                     IntTools_SequenceOfCurves& theSeqOfCurve,
2452                     Standard_Boolean& theTangentFaces,
2453                     Standard_Real& TolReached3d)
2454 {
2455
2456   gp_Pln aPln1 = theS1->Surface().Plane();
2457   gp_Pln aPln2 = theS2->Surface().Plane();
2458
2459   IntAna_QuadQuadGeo aPlnInter(aPln1, aPln2, TolAng, TolTang);
2460
2461   if(!aPlnInter.IsDone()) {
2462     theTangentFaces = Standard_False;
2463     return;
2464   }
2465
2466   IntAna_ResultType aResType = aPlnInter.TypeInter();
2467
2468   if(aResType == IntAna_Same) {
2469     theTangentFaces = Standard_True;
2470     return;
2471   }
2472
2473   theTangentFaces = Standard_False;
2474
2475   if(aResType == IntAna_Empty) {
2476     return;
2477   }
2478
2479   gp_Lin aLin = aPlnInter.Line(1);
2480
2481   ProjLib_Plane aProj;
2482
2483   aProj.Init(aPln1);
2484   aProj.Project(aLin);
2485   gp_Lin2d aLin2d1 = aProj.Line();
2486   //
2487   aProj.Init(aPln2);
2488   aProj.Project(aLin);
2489   gp_Lin2d aLin2d2 = aProj.Line();
2490   //
2491   //classify line2d1 relatively first plane
2492   Standard_Real P11, P12;
2493   Standard_Boolean IsCrossed = ClassifyLin2d(theS1, aLin2d1, TolTang, P11, P12);
2494   if(!IsCrossed) return;
2495   //classify line2d2 relatively second plane
2496   Standard_Real P21, P22;
2497   IsCrossed = ClassifyLin2d(theS2, aLin2d2, TolTang, P21, P22);
2498   if(!IsCrossed) return;
2499
2500   //Analysis of parametric intervals: must have common part
2501
2502   if(P21 >= P12) return;
2503   if(P22 <= P11) return;
2504
2505   Standard_Real pmin, pmax;
2506   pmin = Max(P11, P21);
2507   pmax = Min(P12, P22);
2508
2509   if(pmax - pmin <= TolTang) return;
2510
2511   Handle(Geom_Line) aGLin = new Geom_Line(aLin);
2512
2513   IntTools_Curve aCurve;
2514   Handle(Geom_TrimmedCurve) aGTLin = new Geom_TrimmedCurve(aGLin, pmin, pmax);
2515
2516   aCurve.SetCurve(aGTLin);
2517
2518   if(theApprox1) {
2519     Handle(Geom2d_Line) C2d = new Geom2d_Line(aLin2d1);
2520     aCurve.SetFirstCurve2d(new Geom2d_TrimmedCurve(C2d, pmin, pmax));
2521   }
2522   else {
2523     Handle(Geom2d_Curve) H1;
2524     aCurve.SetFirstCurve2d(H1);
2525   }
2526   if(theApprox2) {
2527     Handle(Geom2d_Line) C2d = new Geom2d_Line(aLin2d2);
2528     aCurve.SetSecondCurve2d(new Geom2d_TrimmedCurve(C2d, pmin, pmax));
2529   }
2530   else {
2531     Handle(Geom2d_Curve) H1;
2532     aCurve.SetFirstCurve2d(H1);
2533   }
2534
2535   theSeqOfCurve.Append(aCurve);
2536   //
2537   // computation of the tolerance reached
2538   Standard_Real anAngle, aDt;
2539   gp_Dir aD1, aD2;
2540   //
2541   aD1 = aPln1.Position().Direction();
2542   aD2 = aPln2.Position().Direction();
2543   anAngle = aD1.Angle(aD2);
2544   //
2545   aDt = IntTools_Tools::ComputeIntRange(TolF1, TolF2, anAngle);
2546   TolReached3d = sqrt(aDt*aDt + TolF1*TolF1);
2547 }
2548
2549 //=======================================================================
2550 //function : ClassifyLin2d
2551 //purpose  :
2552 //=======================================================================
2553 static inline Standard_Boolean INTER(const Standard_Real d1,
2554                                      const Standard_Real d2,
2555                                      const Standard_Real tol)
2556 {
2557   return (d1 > tol && d2 < -tol) ||
2558          (d1 < -tol && d2 > tol) ||
2559          ((d1 <= tol && d1 >= -tol) && (d2 > tol || d2 < -tol)) ||
2560          ((d2 <= tol && d2 >= -tol) && (d1 > tol || d1 < -tol));
2561 }
2562 static inline Standard_Boolean COINC(const Standard_Real d1,
2563                                      const Standard_Real d2,
2564                                      const Standard_Real tol)
2565 {
2566   return (d1 <= tol && d1 >= -tol) && (d2 <= tol && d2 >= -tol);
2567 }
2568 Standard_Boolean ClassifyLin2d(const Handle(GeomAdaptor_HSurface)& theS,
2569                                const gp_Lin2d& theLin2d,
2570                                const Standard_Real theTol,
2571                                Standard_Real& theP1,
2572                                Standard_Real& theP2)
2573
2574 {
2575   Standard_Real xmin, xmax, ymin, ymax, d1, d2, A, B, C;
2576   Standard_Real par[2];
2577   Standard_Integer nbi = 0;
2578
2579   xmin = theS->Surface().FirstUParameter();
2580   xmax = theS->Surface().LastUParameter();
2581   ymin = theS->Surface().FirstVParameter();
2582   ymax = theS->Surface().LastVParameter();
2583
2584   theLin2d.Coefficients(A, B, C);
2585
2586   //xmin, ymin <-> xmin, ymax
2587   d1 = A*xmin + B*ymin + C;
2588   d2 = A*xmin + B*ymax + C;
2589
2590   if(INTER(d1, d2, theTol)) {
2591     //Intersection with boundary
2592     Standard_Real y = -(C + A*xmin)/B;
2593     par[nbi] = ElCLib::Parameter(theLin2d, gp_Pnt2d(xmin, y));
2594     nbi++;
2595   }
2596   else if (COINC(d1, d2, theTol)) {
2597     //Coincidence with boundary
2598     par[0] = ElCLib::Parameter(theLin2d, gp_Pnt2d(xmin, ymin));
2599     par[1] = ElCLib::Parameter(theLin2d, gp_Pnt2d(xmin, ymax));
2600     nbi = 2;
2601   }
2602
2603   if(nbi == 2) {
2604
2605     if(fabs(par[0]-par[1]) > theTol) {
2606       theP1 = Min(par[0], par[1]);
2607       theP2 = Max(par[0], par[1]);
2608       return Standard_True;
2609     }
2610     else return Standard_False;
2611
2612   }
2613
2614   //xmin, ymax <-> xmax, ymax
2615   d1 = d2;
2616   d2 = A*xmax + B*ymax + C;
2617
2618   if(d1 > theTol || d1 < -theTol) {//to avoid checking of
2619                                    //coincidence with the same point
2620     if(INTER(d1, d2, theTol)) {
2621       Standard_Real x = -(C + B*ymax)/A;
2622       par[nbi] = ElCLib::Parameter(theLin2d, gp_Pnt2d(x, ymax));
2623       nbi++;
2624     }
2625     else if (COINC(d1, d2, theTol)) {
2626       par[0] = ElCLib::Parameter(theLin2d, gp_Pnt2d(xmin, ymax));
2627       par[1] = ElCLib::Parameter(theLin2d, gp_Pnt2d(xmax, ymax));
2628       nbi = 2;
2629     }
2630   }
2631
2632   if(nbi == 2) {
2633
2634     if(fabs(par[0]-par[1]) > theTol) {
2635       theP1 = Min(par[0], par[1]);
2636       theP2 = Max(par[0], par[1]);
2637       return Standard_True;
2638     }
2639     else return Standard_False;
2640
2641   }
2642
2643   //xmax, ymax <-> xmax, ymin
2644   d1 = d2;
2645   d2 = A*xmax + B*ymin + C;
2646
2647   if(d1 > theTol || d1 < -theTol) {
2648     if(INTER(d1, d2, theTol)) {
2649       Standard_Real y = -(C + A*xmax)/B;
2650       par[nbi] = ElCLib::Parameter(theLin2d, gp_Pnt2d(xmax, y));
2651       nbi++;
2652     }
2653     else if (COINC(d1, d2, theTol)) {
2654       par[0] = ElCLib::Parameter(theLin2d, gp_Pnt2d(xmax, ymax));
2655       par[1] = ElCLib::Parameter(theLin2d, gp_Pnt2d(xmax, ymin));
2656       nbi = 2;
2657     }
2658   }
2659
2660   if(nbi == 2) {
2661     if(fabs(par[0]-par[1]) > theTol) {
2662       theP1 = Min(par[0], par[1]);
2663       theP2 = Max(par[0], par[1]);
2664       return Standard_True;
2665     }
2666     else return Standard_False;
2667   }
2668
2669   //xmax, ymin <-> xmin, ymin
2670   d1 = d2;
2671   d2 = A*xmin + B*ymin + C;
2672
2673   if(d1 > theTol || d1 < -theTol) {
2674     if(INTER(d1, d2, theTol)) {
2675       Standard_Real x = -(C + B*ymin)/A;
2676       par[nbi] = ElCLib::Parameter(theLin2d, gp_Pnt2d(x, ymin));
2677       nbi++;
2678     }
2679     else if (COINC(d1, d2, theTol)) {
2680       par[0] = ElCLib::Parameter(theLin2d, gp_Pnt2d(xmax, ymin));
2681       par[1] = ElCLib::Parameter(theLin2d, gp_Pnt2d(xmin, ymin));
2682       nbi = 2;
2683     }
2684   }
2685
2686   if(nbi == 2) {
2687     if(fabs(par[0]-par[1]) > theTol) {
2688       theP1 = Min(par[0], par[1]);
2689       theP2 = Max(par[0], par[1]);
2690       return Standard_True;
2691     }
2692     else return Standard_False;
2693   }
2694
2695   return Standard_False;
2696
2697 }
2698 //
2699 //=======================================================================
2700 //function : ApproxParameters
2701 //purpose  :
2702 //=======================================================================
2703 void ApproxParameters(const Handle(GeomAdaptor_HSurface)& aHS1,
2704                       const Handle(GeomAdaptor_HSurface)& aHS2,
2705                       Standard_Integer& iDegMin,
2706                       Standard_Integer& iDegMax,
2707                       Standard_Integer& iNbIter)
2708
2709 {
2710   GeomAbs_SurfaceType aTS1, aTS2;
2711
2712   //
2713   iNbIter=0;
2714   iDegMin=4;
2715   iDegMax=8;
2716   //
2717   aTS1=aHS1->Surface().GetType();
2718   aTS2=aHS2->Surface().GetType();
2719   //
2720   // Cylinder/Torus
2721   if ((aTS1==GeomAbs_Cylinder && aTS2==GeomAbs_Torus) ||
2722       (aTS2==GeomAbs_Cylinder && aTS1==GeomAbs_Torus)) {
2723     Standard_Real aRC, aRT, dR, aPC;
2724     gp_Cylinder aCylinder;
2725     gp_Torus aTorus;
2726     //
2727     aPC=Precision::Confusion();
2728     //
2729     aCylinder=(aTS1==GeomAbs_Cylinder)? aHS1->Surface().Cylinder() : aHS2->Surface().Cylinder();
2730     aTorus=(aTS1==GeomAbs_Torus)? aHS1->Surface().Torus() : aHS2->Surface().Torus();
2731     //
2732     aRC=aCylinder.Radius();
2733     aRT=aTorus.MinorRadius();
2734     dR=aRC-aRT;
2735     if (dR<0.) {
2736       dR=-dR;
2737     }
2738     //
2739     if (dR<aPC) {
2740       iDegMax=6;
2741     }
2742   }
2743   if (aTS1==GeomAbs_Cylinder && aTS2==GeomAbs_Cylinder) {
2744     iNbIter=1;
2745   }
2746 }
2747 //=======================================================================
2748 //function : Tolerances
2749 //purpose  :
2750 //=======================================================================
2751 void Tolerances(const Handle(GeomAdaptor_HSurface)& aHS1,
2752                 const Handle(GeomAdaptor_HSurface)& aHS2,
2753                 Standard_Real& aTolTang)
2754 {
2755   GeomAbs_SurfaceType aTS1, aTS2;
2756   //
2757   aTS1=aHS1->Surface().GetType();
2758   aTS2=aHS2->Surface().GetType();
2759   //
2760   // Cylinder/Torus
2761   if ((aTS1==GeomAbs_Cylinder && aTS2==GeomAbs_Torus) ||
2762       (aTS2==GeomAbs_Cylinder && aTS1==GeomAbs_Torus)) {
2763     Standard_Real aRC, aRT, dR, aPC;
2764     gp_Cylinder aCylinder;
2765     gp_Torus aTorus;
2766     //
2767     aPC=Precision::Confusion();
2768     //
2769     aCylinder=(aTS1==GeomAbs_Cylinder)? aHS1->Surface().Cylinder() : aHS2->Surface().Cylinder();
2770     aTorus=(aTS1==GeomAbs_Torus)? aHS1->Surface().Torus() : aHS2->Surface().Torus();
2771     //
2772     aRC=aCylinder.Radius();
2773     aRT=aTorus.MinorRadius();
2774     dR=aRC-aRT;
2775     if (dR<0.) {
2776       dR=-dR;
2777     }
2778     //
2779     if (dR<aPC) {
2780       aTolTang=0.1*aTolTang;
2781     }
2782   }
2783 }
2784 //=======================================================================
2785 //function : SortTypes
2786 //purpose  :
2787 //=======================================================================
2788 Standard_Boolean SortTypes(const GeomAbs_SurfaceType aType1,
2789                            const GeomAbs_SurfaceType aType2)
2790 {
2791   Standard_Boolean bRet;
2792   Standard_Integer aI1, aI2;
2793   //
2794   bRet=Standard_False;
2795   //
2796   aI1=IndexType(aType1);
2797   aI2=IndexType(aType2);
2798   if (aI1<aI2){
2799     bRet=!bRet;
2800   }
2801   return bRet;
2802 }
2803 //=======================================================================
2804 //function : IndexType
2805 //purpose  :
2806 //=======================================================================
2807 Standard_Integer IndexType(const GeomAbs_SurfaceType aType)
2808 {
2809   Standard_Integer aIndex;
2810   //
2811   aIndex=11;
2812   //
2813   if (aType==GeomAbs_Plane) {
2814     aIndex=0;
2815   }
2816   else if (aType==GeomAbs_Cylinder) {
2817     aIndex=1;
2818   }
2819   else if (aType==GeomAbs_Cone) {
2820     aIndex=2;
2821   }
2822   else if (aType==GeomAbs_Sphere) {
2823     aIndex=3;
2824   }
2825   else if (aType==GeomAbs_Torus) {
2826     aIndex=4;
2827   }
2828   else if (aType==GeomAbs_BezierSurface) {
2829     aIndex=5;
2830   }
2831   else if (aType==GeomAbs_BSplineSurface) {
2832     aIndex=6;
2833   }
2834   else if (aType==GeomAbs_SurfaceOfRevolution) {
2835     aIndex=7;
2836   }
2837   else if (aType==GeomAbs_SurfaceOfExtrusion) {
2838     aIndex=8;
2839   }
2840   else if (aType==GeomAbs_OffsetSurface) {
2841     aIndex=9;
2842   }
2843   else if (aType==GeomAbs_OtherSurface) {
2844     aIndex=10;
2845   }
2846   return aIndex;
2847 }
2848
2849 //=======================================================================
2850 // Function : FindMaxDistance
2851 // purpose :
2852 //=======================================================================
2853 Standard_Real FindMaxDistance(const Handle(Geom_Curve)& theCurve,
2854                               const Standard_Real theFirst,
2855                               const Standard_Real theLast,
2856                               const TopoDS_Face& theFace,
2857                               const Handle(IntTools_Context)& theContext)
2858 {
2859   Standard_Integer aNbS;
2860   Standard_Real aT1, aT2, aDt, aD, aDMax, anEps;
2861   //
2862   aNbS = 11;
2863   aDt = (theLast - theFirst) / aNbS;
2864   aDMax = 0.;
2865   anEps = 1.e-4 * aDt;
2866   //
2867   GeomAPI_ProjectPointOnSurf& aProjPS = theContext->ProjPS(theFace);
2868   aT2 = theFirst;
2869   for (;;) {
2870     aT1 = aT2;
2871     aT2 += aDt;
2872     //
2873     if (aT2 > theLast) {
2874       break;
2875     }
2876     //
2877     aD = FindMaxDistance(theCurve, aT1, aT2, aProjPS, anEps);
2878     if (aD > aDMax) {
2879       aDMax = aD;
2880     }
2881   }
2882   //
2883   return aDMax;
2884 }
2885
2886 //=======================================================================
2887 // Function : FindMaxDistance
2888 // purpose :
2889 //=======================================================================
2890 Standard_Real FindMaxDistance(const Handle(Geom_Curve)& theC,
2891                               const Standard_Real theFirst,
2892                               const Standard_Real theLast,
2893                               GeomAPI_ProjectPointOnSurf& theProjPS,
2894                               const Standard_Real theEps)
2895 {
2896   Standard_Real aA, aB, aCf, aX, aX1, aX2, aF1, aF2, aF;
2897   //
2898   aCf = 0.61803398874989484820458683436564;//(sqrt(5.)-1)/2.;
2899   aA = theFirst;
2900   aB = theLast;
2901   //
2902   aX1=aB - aCf*(aB-aA);
2903   aF1 = MaxDistance(theC, aX1, theProjPS);
2904   aX2 = aA + aCf * (aB - aA);
2905   aF2 = MaxDistance(theC, aX2, theProjPS);
2906
2907   while (Abs(aX1-aX2) > theEps)
2908   {
2909     if (aF1 > aF2) {
2910       aB = aX2;
2911       aX2 = aX1;
2912       aF2 = aF1;
2913       aX1 = aB-aCf*(aB-aA);
2914       aF1 = MaxDistance(theC, aX1, theProjPS);
2915     }
2916     else {
2917       aA = aX1;
2918       aX1 = aX2;
2919       aF1 = aF2;
2920       aX2=aA+aCf*(aB-aA);
2921       aF2 = MaxDistance(theC, aX2, theProjPS);
2922     }
2923   }
2924   //
2925   aX = 0.5 * (aA + aB);
2926   aF = MaxDistance(theC, aX, theProjPS);
2927   //
2928   if (aF1 > aF) {
2929     aF = aF1;
2930   }
2931   //
2932   if (aF2 > aF) {
2933     aF = aF2;
2934   }
2935   //
2936   return aF;
2937 }
2938
2939 //=======================================================================
2940 // Function : MaxDistance
2941 // purpose :
2942 //=======================================================================
2943 Standard_Real MaxDistance(const Handle(Geom_Curve)& theC,
2944                           const Standard_Real aT,
2945                           GeomAPI_ProjectPointOnSurf& theProjPS)
2946 {
2947   Standard_Real aD;
2948   gp_Pnt aP;
2949   //
2950   theC->D0(aT, aP);
2951   theProjPS.Perform(aP);
2952   aD = theProjPS.NbPoints() ? theProjPS.LowerDistance() : 0.;
2953   //
2954   return aD;
2955 }
2956
2957 //=======================================================================
2958 //function : CheckPCurve
2959 //purpose  : Checks if points of the pcurve are out of the face bounds.
2960 //=======================================================================
2961   Standard_Boolean CheckPCurve(const Handle(Geom2d_Curve)& aPC,
2962                                const TopoDS_Face& aFace,
2963                                const Handle(IntTools_Context)& theCtx)
2964 {
2965   const Standard_Integer NPoints = 23;
2966   Standard_Integer i;
2967   Standard_Real umin,umax,vmin,vmax;
2968
2969   theCtx->UVBounds(aFace, umin, umax, vmin, vmax);
2970   Standard_Real tolU = Max ((umax-umin)*0.01, Precision::Confusion());
2971   Standard_Real tolV = Max ((vmax-vmin)*0.01, Precision::Confusion());
2972   Standard_Real fp = aPC->FirstParameter();
2973   Standard_Real lp = aPC->LastParameter();
2974
2975
2976   // adjust domain for periodic surfaces
2977   TopLoc_Location aLoc;
2978   Handle(Geom_Surface) aSurf = BRep_Tool::Surface(aFace, aLoc);
2979   if (aSurf->IsKind(STANDARD_TYPE(Geom_RectangularTrimmedSurface))) {
2980     aSurf = (Handle(Geom_RectangularTrimmedSurface)::DownCast(aSurf))->BasisSurface();
2981   }
2982   gp_Pnt2d pnt = aPC->Value((fp+lp)/2);
2983   Standard_Real u,v;
2984   pnt.Coord(u,v);
2985   //
2986   if (aSurf->IsUPeriodic()) {
2987     Standard_Real aPer = aSurf->UPeriod();
2988     Standard_Integer nshift = (Standard_Integer) ((u-umin)/aPer);
2989     if (u < umin+aPer*nshift) nshift--;
2990     umin += aPer*nshift;
2991     umax += aPer*nshift;
2992   }
2993   if (aSurf->IsVPeriodic()) {
2994     Standard_Real aPer = aSurf->VPeriod();
2995     Standard_Integer nshift = (Standard_Integer) ((v-vmin)/aPer);
2996     if (v < vmin+aPer*nshift) nshift--;
2997     vmin += aPer*nshift;
2998     vmax += aPer*nshift;
2999   }
3000   //
3001   //--------------------------------------------------------
3002   Standard_Boolean bRet;
3003   Standard_Integer j, aNbIntervals;
3004   Standard_Real aT, dT;
3005   gp_Pnt2d aP2D;
3006   //
3007   Geom2dAdaptor_Curve aGAC(aPC);
3008   aNbIntervals=aGAC.NbIntervals(GeomAbs_CN);
3009   //
3010   TColStd_Array1OfReal aTI(1, aNbIntervals+1);
3011   aGAC.Intervals(aTI,GeomAbs_CN);
3012   //
3013   bRet=Standard_False;
3014   //
3015   aT=aGAC.FirstParameter();
3016   for (j=1; j<=aNbIntervals; ++j) {
3017     dT=(aTI(j+1)-aTI(j))/NPoints;
3018     //
3019     for (i=1; i<NPoints; i++) {
3020       aT=aT+dT;
3021       aGAC.D0(aT, aP2D);
3022       aP2D.Coord(u,v);
3023     if (umin-u > tolU || u-umax > tolU ||
3024           vmin-v > tolV || v-vmax > tolV) {
3025         return bRet;
3026   }
3027 }
3028   }
3029   return !bRet;
3030 }
3031 //=======================================================================
3032 //function : CorrectPlaneBoundaries
3033 //purpose  :
3034 //=======================================================================
3035  void CorrectPlaneBoundaries(Standard_Real& aUmin,
3036                              Standard_Real& aUmax,
3037                              Standard_Real& aVmin,
3038                              Standard_Real& aVmax)
3039 {
3040   if (!(Precision::IsInfinite(aUmin) ||
3041         Precision::IsInfinite(aUmax))) {
3042     Standard_Real dU;
3043     //
3044     dU=0.1*(aUmax-aUmin);
3045     aUmin=aUmin-dU;
3046     aUmax=aUmax+dU;
3047   }
3048   if (!(Precision::IsInfinite(aVmin) ||
3049         Precision::IsInfinite(aVmax))) {
3050     Standard_Real dV;
3051     //
3052     dV=0.1*(aVmax-aVmin);
3053     aVmin=aVmin-dV;
3054     aVmax=aVmax+dV;
3055   }
3056 }