src/IntTools/IntTools_EdgeFace.cxx

   1 // Created on: 2001-02-26
   2 // Created by: Peter KURNEV
   3 // Copyright (c) 2001-2012 OPEN CASCADE SAS
   4 //
   5 // The content of this file is subject to the Open CASCADE Technology Public
   6 // License Version 6.5 (the "License"). You may not use the content of this file
   7 // except in compliance with the License. Please obtain a copy of the License
   8 // at http://www.opencascade.org and read it completely before using this file.
   9 //
  10 // The Initial Developer of the Original Code is Open CASCADE S.A.S., having its
  11 // main offices at: 1, place des Freres Montgolfier, 78280 Guyancourt, France.
  12 //
  13 // The Original Code and all software distributed under the License is
  14 // distributed on an "AS IS" basis, without warranty of any kind, and the
  15 // Initial Developer hereby disclaims all such warranties, including without
  16 // limitation, any warranties of merchantability, fitness for a particular
  17 // purpose or non-infringement. Please see the License for the specific terms
  18 // and conditions governing the rights and limitations under the License.
  19
  20
  21 #include <IntTools_EdgeFace.ixx>
  22
  23
  24
  25 #include <IntTools_CArray1OfReal.hxx>
  26 #include <IntTools.hxx>
  27 #include <IntTools_CArray1OfInteger.hxx>
  28 #include <IntTools_Range.hxx>
  29 #include <IntTools_Tools.hxx>
  30 #include <IntTools_Array1OfRange.hxx>
  31 #include <IntTools_QuickSortRange.hxx>
  32 #include <IntTools_CompareRange.hxx>
  33 #include <IntTools_CommonPrt.hxx>
  34 #include <IntTools_Root.hxx>
  35 #include <IntTools_BeanFaceIntersector.hxx>
  36 #include <BOPInt_Context.hxx>
  37
  38 #include <BRep_Tool.hxx>
  39
  40 #include <GeomAdaptor_Surface.hxx>
  41 #include <GeomAdaptor_Curve.hxx>
  42
  43 #include <Geom_Surface.hxx>
  44 #include <Geom_Curve.hxx>
  45
  46 #include <GeomAPI_ProjectPointOnSurf.hxx>
  47
  48 #include <Precision.hxx>
  49
  50 #include <Bnd_Box.hxx>
  51 #include <BndLib_AddSurface.hxx>
  52
  53 #include <gp_Cylinder.hxx>
  54 #include <gp_Ax1.hxx>
  55 #include <gp_Lin.hxx>
  56 #include <gp_Cone.hxx>
  57 #include <gp_Torus.hxx>
  58 #include <gp_Circ.hxx>
  59 #include <gp_Pln.hxx>
  60
  61
  62 #include <Extrema_ExtCS.hxx>
  63 #include <Extrema_POnCurv.hxx>
  64 #include <Extrema_POnSurf.hxx>
  65
  66 // modified by NIZHNY-MKK  Thu Jul 21 11:35:59 2005
  67 #include <IntCurveSurface_HInter.hxx>
  68 #include <GeomAdaptor_HCurve.hxx>
  69 #include <GeomAdaptor_HSurface.hxx>
  70 #include <IntCurveSurface_IntersectionPoint.hxx>
  71
  72 #ifdef WNT
  73 #pragma warning ( disable : 4101 )
  74 #endif
  75
  76 static
  77   Standard_Boolean IsCoplanar (const BRepAdaptor_Curve&  ,
  78                                const BRepAdaptor_Surface& );
  79 static
  80   Standard_Boolean IsRadius (const BRepAdaptor_Curve& aCurve ,
  81                              const BRepAdaptor_Surface& aSurface);
  82 static
  83   Standard_Integer AdaptiveDiscret (const Standard_Integer iDiscret,
  84                                     const BRepAdaptor_Curve& aCurve ,
  85                                     const BRepAdaptor_Surface& aSurface);
  86
  87 //=======================================================================
  88 //function : IntTools_EdgeFace::IntTools_EdgeFace
  89 //purpose  :
  90 //=======================================================================
  91   IntTools_EdgeFace::IntTools_EdgeFace()
  92 {
  93   myTolE=1.e-7;
  94   myTolF=1.e-7;
  95   myDiscret=30;
  96   myEpsT   =1e-12;
  97   myEpsNull=1e-12;
  98   myDeflection=0.01;
  99   myIsDone=Standard_False;
 100   myErrorStatus=1;
 101   myParallel=Standard_False;
 102   myPar1=0.;
 103 }
 104 //=======================================================================
 105 //function : SetContext
 106 //purpose  :
 107 //=======================================================================
 108 void IntTools_EdgeFace::SetContext(const Handle(BOPInt_Context)& theContext)
 109 {
 110   myContext = theContext;
 111 }
 112
 113 //=======================================================================
 114 //function : Context
 115 //purpose  :
 116 //=======================================================================
 117 const Handle(BOPInt_Context)& IntTools_EdgeFace::Context()const
 118 {
 119   return myContext;
 120 }
 121 //=======================================================================
 122 //function : SetEdge
 123 //purpose  :
 124 //=======================================================================
 125   void IntTools_EdgeFace::SetEdge(const TopoDS_Edge& anEdge)
 126 {
 127   myEdge=anEdge;
 128 }
 129
 130 //=======================================================================
 131 //function : SetFace
 132 //purpose  :
 133 //=======================================================================
 134   void IntTools_EdgeFace::SetFace(const TopoDS_Face& aFace)
 135 {
 136   myFace=aFace;
 137 }
 138
 139 //=======================================================================
 140 //function : SetTolE
 141 //purpose  :
 142 //=======================================================================
 143   void IntTools_EdgeFace::SetTolE(const Standard_Real aTol)
 144 {
 145   myTolE=aTol;
 146 }
 147 //=======================================================================
 148 //function : SetTolF
 149 //purpose  :
 150 //=======================================================================
 151   void IntTools_EdgeFace::SetTolF(const Standard_Real aTol)
 152 {
 153   myTolF=aTol;
 154 }
 155
 156 //=======================================================================
 157 //function : SetDiscretize
 158 //purpose  :
 159 //=======================================================================
 160   void IntTools_EdgeFace::SetDiscretize(const Standard_Integer aDiscret)
 161 {
 162   myDiscret=aDiscret;
 163 }
 164 //=======================================================================
 165 //function : SetDeflection
 166 //purpose  :
 167 //=======================================================================
 168   void IntTools_EdgeFace::SetDeflection(const Standard_Real aDefl)
 169 {
 170   myDeflection=aDefl;
 171 }
 172 //=======================================================================
 173 //function : SetEpsilonT
 174 //purpose  :
 175 //=======================================================================
 176   void IntTools_EdgeFace::SetEpsilonT(const Standard_Real anEpsT)
 177 {
 178   myEpsT=anEpsT;
 179 }
 180 //=======================================================================
 181 //function : SetEpsilonNull
 182 //purpose  :
 183 //=======================================================================
 184   void IntTools_EdgeFace::SetEpsilonNull(const Standard_Real anEpsNull)
 185 {
 186   myEpsNull=anEpsNull;
 187 }
 188
 189 //=======================================================================
 190 //function : SetRange
 191 //purpose  :
 192 //=======================================================================
 193   void IntTools_EdgeFace::SetRange(const Standard_Real aFirst,
 194                                    const Standard_Real aLast)
 195 {
 196   myRange.SetFirst (aFirst);
 197   myRange.SetLast  (aLast);
 198 }
 199
 200 //=======================================================================
 201 //function : SetRange
 202 //purpose  :
 203 //=======================================================================
 204   void IntTools_EdgeFace::SetRange(const IntTools_Range& aRange)
 205 {
 206   myRange.SetFirst (aRange.First());
 207   myRange.SetLast  (aRange.Last());
 208 }
 209 //=======================================================================
 210 //function : IsDone
 211 //purpose  :
 212 //=======================================================================
 213 Standard_Boolean IntTools_EdgeFace::IsDone()const
 214 {
 215   return myIsDone;
 216 }
 217 //=======================================================================
 218 //function : ErrorStatus
 219 //purpose  :
 220 //=======================================================================
 221 Standard_Integer IntTools_EdgeFace::ErrorStatus()const
 222 {
 223   return myErrorStatus;
 224 }
 225 //=======================================================================
 226 //function : CommonParts
 227 //purpose  :
 228 //=======================================================================
 229 const IntTools_SequenceOfCommonPrts& IntTools_EdgeFace::CommonParts() const
 230 {
 231   return mySeqOfCommonPrts;
 232 }
 233 //=======================================================================
 234 //function : Range
 235 //purpose  :
 236 //=======================================================================
 237 const IntTools_Range&  IntTools_EdgeFace::Range() const
 238 {
 239   return myRange;
 240 }
 241
 242 //=======================================================================
 243 //function : CheckData
 244 //purpose  :
 245 //=======================================================================
 246 void IntTools_EdgeFace::CheckData()
 247 {
 248   if (BRep_Tool::Degenerated(myEdge)) {
 249     myErrorStatus=2;
 250   }
 251   if (!BRep_Tool::IsGeometric(myEdge)) {
 252      myErrorStatus=3;
 253   }
 254 }
 255 //=======================================================================
 256 //function : Prepare
 257 //purpose  :
 258 //=======================================================================
 259   void IntTools_EdgeFace::Prepare()
 260 {
 261   Standard_Integer pri;
 262   IntTools_CArray1OfReal aPars;
 263
 264   //
 265   // 1.Prepare Curve's data and Surface's data
 266   myC.Initialize(myEdge);
 267   GeomAbs_CurveType aCurveType;
 268   aCurveType=myC.GetType();
 269   //
 270   // 2.Prepare myCriteria
 271   if (aCurveType==GeomAbs_BSplineCurve||
 272         aCurveType==GeomAbs_BezierCurve) {
 273     myCriteria=1.5*myTolE+myTolF;
 274   }
 275   else {
 276     myCriteria=myTolE+myTolF;
 277   }
 278   // 2.a myTmin, myTmax
 279   myTmin=myRange.First();
 280   myTmax=myRange.Last();
 281   // 2.b myFClass2d
 282   myS.Initialize (myFace,Standard_True);
 283   myFClass2d.Init(myFace, 1.e-6);
 284   //
 285   // 2.c Prepare adaptive myDiscret
 286   myDiscret=AdaptiveDiscret(myDiscret, myC, myS);
 287   //
 288   //
 289   // 3.Prepare myPars
 290   pri = IntTools::PrepareArgs(myC, myTmax, myTmin, myDiscret, myDeflection, aPars);
 291   if (pri) {
 292     myErrorStatus=6;
 293     return;
 294   }
 295   // 4.
 296   //ProjectableRanges
 297   Standard_Integer i, iProj, aNb, aNbProj, ind0, ind1;
 298   Standard_Real t0, t1, tRoot;
 299
 300   //
 301   // Table of Projection's function values
 302   aNb=aPars.Length();
 303   IntTools_CArray1OfInteger anArrProjectability;
 304   anArrProjectability.Resize(aNb);
 305
 306   for (iProj=0, i=0; i<aNb; i++) {
 307     t0=aPars(i);
 308     aNbProj=IsProjectable (t0);
 309
 310     anArrProjectability(i)=0;
 311     if (aNbProj) {
 312       anArrProjectability(i)=1;
 313       iProj++;
 314     }
 315   }
 316   //
 317   // Checking
 318   if (!iProj ) {
 319     myErrorStatus=7;
 320     return;
 321   }
 322
 323   //
 324   // Projectable Ranges
 325   IntTools_Range aRange;
 326
 327   ind0=anArrProjectability(0);
 328   if (ind0) {
 329     t0=aPars(0);
 330     aRange.SetFirst(t0);
 331   }
 332
 333   for(i=1; i<aNb; i++) {
 334     ind1=anArrProjectability(i);
 335     t0=aPars(i-1);
 336     t1=aPars(i);
 337
 338     if (i==(aNb-1)) {
 339       if (ind1 && ind0) {
 340         aRange.SetLast(t1);
 341         myProjectableRanges.Append(aRange);
 342       }
 343       if (ind1 && !ind0) {
 344         FindProjectableRoot(t0, t1, ind0, ind1, tRoot);
 345         aRange.SetFirst(tRoot);
 346         aRange.SetLast(t1);
 347         myProjectableRanges.Append(aRange);
 348       }
 349       //
 350       if (ind0 && !ind1) {
 351         FindProjectableRoot(t0, t1, ind0, ind1, tRoot);
 352         aRange.SetLast(tRoot);
 353         myProjectableRanges.Append(aRange);
 354       }
 355       //
 356       break;
 357     }
 358
 359     if (ind0 != ind1) {
 360       FindProjectableRoot(t0, t1, ind0, ind1, tRoot);
 361
 362       if (ind0 && !ind1) {
 363         aRange.SetLast(tRoot);
 364         myProjectableRanges.Append(aRange);
 365       }
 366       else {
 367         aRange.SetFirst(tRoot);
 368       }
 369     } // if (ind0 != ind1)
 370     ind0=ind1;
 371   } // for(i=1; i<aNb; i++) {
 372 }
 373
 374 //=======================================================================
 375 //function : FindProjectableRoot
 376 //purpose  :
 377 //=======================================================================
 378   void IntTools_EdgeFace::FindProjectableRoot (const Standard_Real tt1,
 379                                                const Standard_Real tt2,
 380                                                const Standard_Integer ff1,
 381                                                const Standard_Integer ff2,
 382                                                Standard_Real& tRoot)
 383 {
 384   Standard_Real tm, t1, t2, aEpsT;
 385   Standard_Integer anIsProj1, anIsProj2, anIsProjm;
 386   aEpsT=0.5*myEpsT;
 387   //
 388   // Root is inside [tt1, tt2]
 389   t1=tt1;
 390   t2=tt2;
 391   anIsProj1=ff1;
 392   anIsProj2=ff2;
 393
 394   while (1) {
 395     if (fabs(t1-t2) < aEpsT) {
 396       tRoot=(anIsProj1) ? t1 : t2;
 397       return;
 398     }
 399     tm=.5*(t1+t2);
 400     anIsProjm=IsProjectable(tm);
 401
 402     if (anIsProjm != anIsProj1) {
 403       t2=tm;
 404       anIsProj2=anIsProjm;
 405     }
 406     else {
 407       t1=tm;
 408       anIsProj1=anIsProjm;
 409     }
 410   }
 411 }
 412 //=======================================================================
 413 //function : IsProjectable
 414 //purpose  :
 415 //=======================================================================
 416   Standard_Boolean IntTools_EdgeFace::IsProjectable(const Standard_Real t) const
 417 {
 418   Standard_Boolean bFlag;
 419   Standard_Real Umin, Usup, Vmin, Vsup;
 420
 421   Umin=myS.FirstUParameter();
 422   Usup=myS.LastUParameter();
 423   Vmin=myS.FirstVParameter();
 424   Vsup=myS.LastVParameter ();
 425
 426   gp_Pnt P;
 427   myC.D0(t, P);
 428   GeomAPI_ProjectPointOnSurf aProjector;
 429   //
 430   Standard_Real ULD, VLD;
 431
 432   GeomAPI_ProjectPointOnSurf& aLocProj = myContext->ProjPS(myFace);
 433   aLocProj.Perform(P);
 434   bFlag = aLocProj.IsDone();
 435
 436   if(bFlag) {
 437     aLocProj.LowerDistanceParameters(ULD, VLD);
 438   }
 439   //
 440
 441   if (bFlag) {
 442     bFlag=Standard_False;
 443
 444     //
 445     TopAbs_State aState;
 446     gp_Pnt2d aP2d(ULD, VLD);
 447
 448     aState = myContext->FClass2d(myFace).Perform(aP2d);
 449     //
 450
 451     if (aState==TopAbs_IN || aState==TopAbs_ON) {
 452       bFlag=Standard_True;
 453     }
 454   }
 455   return bFlag;
 456 }
 457
 458 //=======================================================================
 459 //function : DistanceFunction
 460 //purpose  :
 461 //=======================================================================
 462   Standard_Real IntTools_EdgeFace::DistanceFunction(const Standard_Real t)
 463 {
 464   Standard_Real Umin, Usup, Vmin, Vsup, aD;
 465   //
 466   gp_Pnt P;
 467   myC.D0(t, P);
 468   //
 469   Standard_Boolean bIsEqDistance;
 470
 471   bIsEqDistance= IntTools_EdgeFace::IsEqDistance(P, myS, 1.e-7, aD);
 472   if (bIsEqDistance) {
 473     aD=aD-myCriteria;
 474     return aD;
 475   }
 476
 477   Umin=myS.FirstUParameter();
 478   Usup=myS.LastUParameter();
 479   Vmin=myS.FirstVParameter();
 480   Vsup=myS.LastVParameter ();
 481
 482   //
 483   Standard_Boolean bFlag = Standard_False;
 484
 485   GeomAPI_ProjectPointOnSurf& aLocProj = myContext->ProjPS(myFace);
 486   aLocProj.Perform(P);
 487   bFlag = aLocProj.IsDone();
 488
 489   if(bFlag) {
 490     aD = aLocProj.LowerDistance();
 491   }
 492   //
 493
 494   if (!bFlag) {
 495     myErrorStatus=11;
 496     return 99.;
 497   }
 498
 499   //
 500   //   aD=aProjector.LowerDistance();
 501   //
 502   aD=aD-myCriteria;
 503   return aD;
 504 }
 505 //
 506 //=======================================================================
 507 //function : IsEqDistance
 508 //purpose  :
 509 //=======================================================================
 510   Standard_Boolean IntTools_EdgeFace::IsEqDistance(const gp_Pnt& aP,
 511                                                    const BRepAdaptor_Surface& aBAS,
 512                                                    const Standard_Real aTol,
 513                                                    Standard_Real& aD)
 514 {
 515   Standard_Boolean bRetFlag=Standard_True;
 516
 517   GeomAbs_SurfaceType aSurfType=aBAS.GetType();
 518
 519   if (aSurfType==GeomAbs_Cylinder) {
 520     gp_Cylinder aCyl=aBAS.Cylinder();
 521     const gp_Ax1& anAx1  =aCyl.Axis();
 522     gp_Lin aLinAxis(anAx1);
 523     Standard_Real aDC, aRadius=aCyl.Radius();
 524     aDC=aLinAxis.Distance(aP);
 525     if (aDC < aTol) {
 526       aD=aRadius;
 527       return bRetFlag;
 528     }
 529   }
 530
 531   if (aSurfType==GeomAbs_Cone) {
 532     gp_Cone aCone=aBAS.Cone();
 533     const gp_Ax1& anAx1  =aCone.Axis();
 534     gp_Lin aLinAxis(anAx1);
 535     Standard_Real aDC, aRadius, aDS, aSemiAngle;
 536     aDC=aLinAxis.Distance(aP);
 537     if (aDC < aTol) {
 538       gp_Pnt anApex=aCone.Apex();
 539       aSemiAngle=aCone.SemiAngle();
 540       aDS=aP.Distance(anApex);
 541
 542       aRadius=aDS*tan(aSemiAngle);
 543       aD=aRadius;
 544       return bRetFlag;
 545     }
 546   }
 547
 548   if (aSurfType==GeomAbs_Torus) {
 549     Standard_Real aMajorRadius, aMinorRadius, aDC;
 550
 551     gp_Torus aTorus=aBAS.Torus();
 552     gp_Pnt aPLoc=aTorus.Location();
 553     aMajorRadius=aTorus.MajorRadius();
 554
 555     aDC=fabs(aPLoc.Distance(aP)-aMajorRadius);
 556     if (aDC < aTol) {
 557       aMinorRadius=aTorus.MinorRadius();
 558       aD=aMinorRadius;
 559       return bRetFlag;
 560     }
 561   }
 562   return !bRetFlag;
 563 }
 564 //
 565 //=======================================================================
 566 //function : PrepareArgsFuncArrays
 567 //purpose  : Obtain
 568 //           myFuncArray and myArgsArray for the interval [ta, tb]
 569 //=======================================================================
 570   void IntTools_EdgeFace::PrepareArgsFuncArrays(const Standard_Real ta,
 571                                                 const Standard_Real tb)
 572 {
 573   IntTools_CArray1OfReal anArgs, aFunc;
 574   Standard_Integer i, aNb, pri;
 575   Standard_Real t, f, f1;
 576   //
 577   // Prepare values of arguments for the interval [ta, tb]
 578   pri=IntTools::PrepareArgs (myC, tb, ta, myDiscret, myDeflection, anArgs);
 579
 580   if (pri) {
 581     myErrorStatus=8;
 582     return;
 583   }
 584   //...
 585   aNb=anArgs.Length();
 586
 587   if (!aNb){
 588     myErrorStatus=9;
 589     return;
 590   }
 591   //
 592   // Prepare values of functions for the interval [ta, tb]
 593   aFunc.Resize(aNb);
 594   for (i=0; i<aNb; i++) {
 595     t=anArgs(i);
 596     f1=DistanceFunction(t);
 597     f=f1+myCriteria;
 598
 599     if (myErrorStatus==11)
 600       return;
 601
 602     if (f1 < myEpsNull) {
 603       f=0.;
 604     }
 605     aFunc(i)=f;
 606   }
 607   //
 608   // Add points where the derivative = 0
 609   AddDerivativePoints(anArgs, aFunc);
 610
 611 }
 612
 613
 614 //=======================================================================
 615 //function : AddDerivativePoints
 616 //purpose  :
 617 //=======================================================================
 618   void IntTools_EdgeFace::AddDerivativePoints(const IntTools_CArray1OfReal& t,
 619                                               const IntTools_CArray1OfReal& f)
 620 {
 621   Standard_Integer i, j, n, k, nn=100;
 622   Standard_Real fr, tr, tr1, dEpsNull=10.*myEpsNull;
 623   IntTools_CArray1OfReal fd;
 624   TColStd_SequenceOfReal aTSeq, aFSeq;
 625
 626   n=t.Length();
 627   fd.Resize(n+1);
 628   //
 629   // Table of derivatives
 630   Standard_Real dfx, tx, tx1, fx, fx1, dt=1.e-6;
 631   // Left limit
 632   tx=t(0);
 633   tx1=tx+dt;
 634   fx=f(0);
 635   fx1=DistanceFunction(tx1);
 636   fx1=fx1+myCriteria;
 637   if (fx1 < myEpsNull) {
 638     fx1=0.;
 639   }
 640   dfx=(fx1-fx)/dt;
 641   fd(0)=dfx;
 642
 643   if (fabs(fd(0)) < dEpsNull){
 644     fd(0)=0.;
 645   }
 646
 647
 648   k=n-1;
 649   for (i=1; i<k; i++) {
 650     Standard_Real ti, ti1;
 651     ti=t(i);
 652     ti1=t(i-1);
 653     fd(i)=.5*(f(i+1)-f(i-1))/(t(i)-t(i-1));
 654     if (fabs(fd(i)) < dEpsNull){
 655       fd(i)=0.;
 656     }
 657   }
 658   // Right limit
 659   tx=t(n-1);
 660   tx1=tx-dt;
 661   fx=f(n-1);
 662   fx1=DistanceFunction(tx1);
 663   fx1=fx1+myCriteria;
 664   if (fx1 < myEpsNull) {
 665     fx1=0.;
 666   }
 667   dfx=(fx-fx1)/dt;
 668   fd(n-1)=dfx;
 669
 670   if (fabs(fd(n-1)) < dEpsNull){
 671     fd(n-1)=0.;
 672   }
 673   //
 674   // Finding the range where the derivatives have different signs
 675   // for neighbouring points
 676   for (i=1; i<n; i++) {
 677     Standard_Real fd1, fd2, t1, t2;
 678     t1 =t(i-1);
 679     t2 =t(i);
 680     fd1=fd(i-1);
 681     fd2=fd(i);
 682
 683     if (fd1*fd2 < 0.) {
 684       if (fabs(fd1) < myEpsNull) {
 685         tr=t1;
 686         fr=DistanceFunction(tr);//fd1;
 687       }
 688       else if (fabs(fd2) < myEpsNull) {
 689         tr=t2;
 690         fr=DistanceFunction(tr);
 691       }
 692       else {
 693         tr=FindSimpleRoot(2, t1, t2, fd1);
 694         fr=DistanceFunction(tr);
 695       }
 696
 697       aTSeq.Append(tr);
 698       aFSeq.Append(fr);
 699     }
 700   } // end of for (i=1; i<n; i++)
 701   //
 702   // remove identical t, f
 703   nn=aTSeq.Length();
 704   if (nn) {
 705     for (i=1; i<=aTSeq.Length(); i++) {
 706       tr=aTSeq(i);
 707       for (j=0; j<n; j++) {
 708         tr1=t(j);
 709         if (fabs (tr1-tr) < myEpsT) {
 710           aTSeq.Remove(i);
 711           aFSeq.Remove(i);
 712         }
 713       }
 714     }
 715     nn=aTSeq.Length();
 716   }
 717   //
 718   // sorting args and funcs in increasing order
 719   if (nn) {
 720     k=nn+n;
 721     IntTools_Array1OfRange anArray1OfRange(1, k);
 722     for (i=1; i<=n; i++) {
 723       anArray1OfRange(i).SetFirst(t(i-1));
 724       anArray1OfRange(i).SetLast (f(i-1));
 725     }
 726     for (i=1; i<=nn; i++) {
 727       anArray1OfRange(n+i).SetFirst(aTSeq(i));
 728       anArray1OfRange(n+i).SetLast (aFSeq(i));
 729     }
 730
 731     IntTools_QuickSortRange aQuickSortRange;
 732     IntTools_CompareRange aComparator;
 733     aQuickSortRange.Sort (anArray1OfRange, aComparator);
 734
 735     // filling the  output arrays
 736     myArgsArray.Resize(k);
 737     myFuncArray.Resize(k);
 738     for (i=1; i<=k; i++) {
 739       myArgsArray(i-1)=anArray1OfRange(i).First();
 740       myFuncArray(i-1)=anArray1OfRange(i).Last ();
 741     }
 742   }
 743
 744   else { // nn=0
 745     myArgsArray.Resize(n);
 746     myFuncArray.Resize(n);
 747     for (i=0; i<n; i++) {
 748       myArgsArray(i)=t(i);
 749       myFuncArray(i)=f(i);
 750     }
 751   }
 752 }
 753
 754 //=======================================================================
 755 //function : DerivativeFunction
 756 //purpose  :
 757 //=======================================================================
 758   Standard_Real IntTools_EdgeFace::DerivativeFunction(const Standard_Real t2)
 759 {
 760   Standard_Real t1, t3, aD1, aD2, aD3;
 761   Standard_Real dt=1.e-9;
 762   t1=t2-dt;
 763   aD1=DistanceFunction(t1);
 764   t3=t2+dt;
 765   aD3=DistanceFunction(t3);
 766
 767   aD2=.5*(aD3-aD1)/dt;
 768   return aD2;
 769 }
 770
 771 //=======================================================================
 772 //function : FindSimpleRoot
 773 //purpose  : [private]
 774 //=======================================================================
 775   Standard_Real IntTools_EdgeFace::FindSimpleRoot (const Standard_Integer IP,
 776                                                    const Standard_Real tA,
 777                                                    const Standard_Real tB,
 778                                                    const Standard_Real fA)
 779 {
 780   Standard_Real r, a, b, y, x0, s;
 781
 782   a=tA; b=tB; r=fA;
 783
 784   while (1) {
 785     x0=.5*(a+b);
 786
 787     if (IP==1)
 788       y=DistanceFunction(x0);
 789     else
 790       y=DerivativeFunction(x0);
 791
 792     if (fabs(b-a) < myEpsT || y==0.) {
 793       return x0;
 794     }
 795
 796
 797     s=y*r;
 798
 799     if (s<0.) {
 800       b=x0;
 801       continue;
 802     }
 803
 804     if (s>0.) {
 805       a=x0; r=y;
 806     }
 807   }
 808 }
 809 //=======================================================================
 810 //function : FindGoldRoot
 811 //purpose  : [private]
 812 //=======================================================================
 813   Standard_Real IntTools_EdgeFace::FindGoldRoot (const Standard_Real tA,
 814                                                  const Standard_Real tB,
 815                                                  const Standard_Real coeff)
 816 {
 817   Standard_Real gs=0.61803399;
 818   Standard_Real a, b, xp, xl, yp, yl;
 819
 820   a=tA;  b=tB;
 821
 822   xp=a+(b-a)*gs;
 823   xl=b-(b-a)*gs;
 824   yp=coeff*DistanceFunction(xp);
 825   yl=coeff*DistanceFunction(xl);
 826
 827
 828   while (1) {
 829
 830     if (fabs(b-a) < myEpsT) {
 831       return .5*(b+a);
 832     }
 833
 834     if (yp < yl) {
 835       a=xl;
 836       xl=xp;
 837       xp=a+(b-a)*gs;
 838       yp=coeff*DistanceFunction(xp);
 839     }
 840
 841     else {
 842       b=xp;
 843       xp=xl;
 844       yp=yl;
 845       xl=b-(b-a)*gs;
 846       yl=coeff*DistanceFunction(xl);
 847     }
 848   }
 849 }
 850
 851 //=======================================================================
 852 //function : MakeType
 853 //purpose  :
 854 //=======================================================================
 855   Standard_Integer IntTools_EdgeFace::MakeType(IntTools_CommonPrt&  aCommonPrt)
 856 {
 857   Standard_Real  af1, al1;
 858   Standard_Real dt, df1, df2, tm;
 859   Standard_Boolean bAllNullFlag;
 860   //
 861   bAllNullFlag=aCommonPrt.AllNullFlag();
 862   if (bAllNullFlag) {
 863     aCommonPrt.SetType(TopAbs_EDGE);
 864     return 0;
 865   }
 866   //
 867   aCommonPrt.Range1(af1, al1);
 868
 869   {
 870     gp_Pnt aPF, aPL;
 871     myC.D0(af1, aPF);
 872     myC.D0(al1, aPL);
 873     df1=aPF.Distance(aPL);
 874     Standard_Boolean isWholeRange = Standard_False;
 875
 876     if((Abs(af1 - myRange.First()) < myC.Resolution(myCriteria)) &&
 877        (Abs(al1 - myRange.Last()) < myC.Resolution(myCriteria)))
 878       isWholeRange = Standard_True;
 879
 880
 881     if ((df1 > myCriteria * 2.) && isWholeRange) {
 882       aCommonPrt.SetType(TopAbs_EDGE);
 883     }
 884     else {
 885       if(isWholeRange) {
 886         tm = (af1 + al1) * 0.5;
 887
 888         if(aPF.Distance(myC.Value(tm)) > myCriteria * 2.) {
 889           aCommonPrt.SetType(TopAbs_EDGE);
 890           return 0;
 891         }
 892       }
 893
 894       if(!CheckTouch(aCommonPrt, tm)) {
 895         tm = (af1 + al1) * 0.5;
 896       }
 897       aCommonPrt.SetType(TopAbs_VERTEX);
 898       aCommonPrt.SetVertexParameter1(tm);
 899       aCommonPrt.SetRange1 (af1, al1);
 900     }
 901          return 0;
 902   }
 903   //
 904   /*
 905   dt=al1-af1;
 906   if (dt<1.e-5) {
 907     gp_Pnt aPF, aPL;
 908     myC.D0(af1, aPF);
 909     myC.D0(al1, aPL);
 910     df1=aPF.Distance(aPL);
 911     if (df1<myCriteria) {
 912       //
 913       tm=.5*(af1+al1);
 914       aCommonPrt.SetType(TopAbs_VERTEX);
 915       aCommonPrt.SetVertexParameter1(tm);
 916       aCommonPrt.SetRange1 (af1, al1);
 917       return 0;
 918     }
 919   }
 920   //
 921   IsIntersection (af1, al1);
 922   //
 923   if (!myParallel) {
 924     aCommonPrt.SetType(TopAbs_VERTEX);
 925     aCommonPrt.SetVertexParameter1(myPar1);
 926     aCommonPrt.SetRange1 (af1, al1);
 927   }
 928   else {
 929     dt=al1-af1;
 930     if (dt<1.e-5) {
 931       df1=DistanceFunction(af1);
 932       df2=DistanceFunction(al1);
 933       tm=(df1 < df2) ? af1 : al1;
 934       aCommonPrt.SetType(TopAbs_VERTEX);
 935       aCommonPrt.SetVertexParameter1(tm);
 936       aCommonPrt.SetRange1 (af1, al1);
 937     }
 938
 939     else {
 940       aCommonPrt.SetType(TopAbs_EDGE);
 941     }
 942   }
 943   */
 944   return 0;
 945 }
 946
 947
 948 //=======================================================================
 949 //function : IsIntersection
 950 //purpose  :
 951 //=======================================================================
 952   void IntTools_EdgeFace::IsIntersection (const Standard_Real ta,
 953                                           const Standard_Real tb)
 954 {
 955   IntTools_CArray1OfReal anArgs, aFunc;
 956   Standard_Integer i, aNb, pri, aCnt=0;
 957   //
 958   Standard_Integer aCntIncreasing=1, aCntDecreasing=1;
 959   Standard_Real t, f, f1;
 960   //
 961   // Prepare values of arguments for the interval [ta, tb]
 962   pri=IntTools::PrepareArgs (myC, tb, ta, myDiscret, myDeflection, anArgs);
 963   aNb=anArgs.Length();
 964
 965   aFunc.Resize(aNb);
 966   for (i=0; i<aNb; i++) {
 967     t=anArgs(i);
 968
 969     f1=DistanceFunction(t);
 970     f=f1+myCriteria;
 971
 972     if (fabs(f1) < myEpsNull) {
 973       aCnt++;
 974       f=0.;
 975     }
 976     aFunc(i)=f;
 977     //
 978     if (i) {
 979       if (aFunc(i)>aFunc(i-1)) {
 980         aCntIncreasing++;
 981       }
 982       if (aFunc(i)<aFunc(i-1)) {
 983         aCntDecreasing++;
 984       }
 985     }
 986     //
 987   }
 988
 989   if (aCnt==aNb) {
 990     myParallel=Standard_True;
 991     return;
 992   }
 993
 994   FindDerivativeRoot(anArgs, aFunc);
 995
 996   //
 997   if (myParallel) {
 998     if (!(myC.GetType()==GeomAbs_Line
 999           &&
1000           myS.GetType()==GeomAbs_Cylinder)) {
1001       if (aCntDecreasing==aNb) {
1002         myPar1=anArgs(aNb-1);
1003         myParallel=Standard_False;
1004       }
1005       if (aCntIncreasing==aNb) {
1006         myPar1=anArgs(0);
1007         myParallel=Standard_False;
1008       }
1009     }
1010   }
1011   //
1012   return ;
1013 }
1014
1015 //=======================================================================
1016 //function : FindDerivativeRoot
1017 //purpose  :
1018 //=======================================================================
1019   void IntTools_EdgeFace::FindDerivativeRoot(const IntTools_CArray1OfReal& t,
1020                                              const IntTools_CArray1OfReal& f)
1021 {
1022   Standard_Integer i, n, k;
1023   Standard_Real fr, tr;
1024   IntTools_CArray1OfReal fd;
1025   TColStd_SequenceOfReal aTSeq, aFSeq;
1026
1027   myPar1=0.;
1028   myParallel=Standard_True;
1029
1030   n=t.Length();
1031   fd.Resize(n+1);
1032   //
1033   // Table of derivatives
1034   fd(0)=(f(1)-f(0))/(t(1)-t(0));
1035   if (fabs(fd(0)) < myEpsNull) {
1036     fd(0)=0.;
1037   }
1038
1039   k=n-1;
1040   for (i=1; i<k; i++) {
1041     fd(i)=.5*(f(i+1)-f(i-1))/(t(i)-t(i-1));
1042     if (fabs(fd(i)) < myEpsNull) {
1043       fd(i)=0.;
1044     }
1045   }
1046
1047   fd(n-1)=(f(n-1)-f(n-2))/(t(n-1)-t(n-2));
1048   if (fabs(fd(n-1)) < myEpsNull) {
1049     fd(n-1)=0.;
1050   }
1051   //
1052   // Finding the range where the derivatives have different signs
1053   // for neighbouring points
1054   for (i=1; i<n; i++) {
1055     Standard_Real fd1, fd2, t1, t2, fabsfd1, fabsfd2;
1056     Standard_Boolean bF1, bF2;
1057     t1 =t(i-1);
1058     t2 =t(i);
1059     fd1=fd(i-1);
1060     fd2=fd(i);
1061
1062     fabsfd1=fabs(fd1);
1063     bF1=fabsfd1 < myEpsNull;
1064
1065     fabsfd2=fabs(fd2);
1066     bF2=fabsfd2 < myEpsNull;
1067     //
1068     if (fd1*fd2 < 0.) {
1069       tr=FindSimpleRoot(2, t1, t2, fd1);
1070       fr=DistanceFunction(tr);
1071       myPar1=tr;
1072       myParallel=Standard_False;
1073       break;
1074     }
1075
1076     if (!bF1 && bF2) {
1077       tr=t2;
1078       fr=fd2;
1079       myPar1=tr;
1080       myParallel=Standard_False;
1081       break;
1082     }
1083
1084     if (bF1 && !bF2) {
1085       tr=t1;
1086       fr=fd1;
1087       myPar1=tr;
1088       myParallel=Standard_False;
1089       break;
1090     }
1091
1092   }
1093 }
1094 //=======================================================================
1095 //function : RemoveIdenticalRoots
1096 //purpose  :
1097 //=======================================================================
1098   void IntTools_EdgeFace::RemoveIdenticalRoots()
1099 {
1100   Standard_Integer aNbRoots, j, k;
1101
1102   aNbRoots=mySequenceOfRoots.Length();
1103   for (j=1; j<=aNbRoots; j++) {
1104     const IntTools_Root& aRj=mySequenceOfRoots(j);
1105     for (k=j+1; k<=aNbRoots; k++) {
1106       const IntTools_Root& aRk=mySequenceOfRoots(k);
1107
1108       Standard_Real aTj, aTk, aDistance;
1109       gp_Pnt aPj, aPk;
1110
1111       aTj=aRj.Root();
1112       aTk=aRk.Root();
1113
1114       myC.D0(aTj, aPj);
1115       myC.D0(aTk, aPk);
1116
1117       aDistance=aPj.Distance(aPk);
1118       if (aDistance < myCriteria) {
1119         mySequenceOfRoots.Remove(k);
1120         aNbRoots=mySequenceOfRoots.Length();
1121       }
1122     }
1123   }
1124 }
1125
1126 //=======================================================================
1127 //function : CheckTouch
1128 //purpose  :
1129 //=======================================================================
1130   Standard_Boolean IntTools_EdgeFace::CheckTouch(const IntTools_CommonPrt& aCP,
1131                                                  Standard_Real&            aTx)
1132 {
1133   Standard_Real aTF, aTL, Tol, U1f, U1l, V1f, V1l, af, al,aDist2, aMinDist2;
1134   Standard_Boolean theflag=Standard_False;
1135   Standard_Integer aNbExt, i, iLower ;
1136
1137   aCP.Range1(aTF, aTL);
1138
1139   //
1140   Standard_Real aCR;
1141   aCR=myC.Resolution(myCriteria);
1142   if((Abs(aTF - myRange.First()) < aCR) &&
1143      (Abs(aTL - myRange.Last())  < aCR)) {
1144     return theflag; // EDGE
1145   }
1146   //
1147
1148   Tol = Precision::PConfusion();
1149
1150   const Handle(Geom_Curve)&  Curve   =BRep_Tool::Curve  (myC.Edge(), af, al);
1151   const Handle(Geom_Surface)& Surface=BRep_Tool::Surface(myS.Face());
1152   //   Surface->Bounds(U1f,U1l,V1f,V1l);
1153   U1f = myS.FirstUParameter();
1154   U1l = myS.LastUParameter();
1155   V1f = myS.FirstVParameter();
1156   V1l = myS.LastVParameter();
1157
1158   GeomAdaptor_Curve   TheCurve   (Curve,aTF, aTL);
1159   GeomAdaptor_Surface TheSurface (Surface, U1f, U1l, V1f, V1l);
1160
1161   Extrema_ExtCS anExtrema (TheCurve, TheSurface, Tol, Tol);
1162
1163   aDist2 = 1.e100;
1164
1165   if(anExtrema.IsDone()) {
1166     aMinDist2 = aDist2;
1167
1168     if(!anExtrema.IsParallel()) {
1169       aNbExt=anExtrema.NbExt();
1170
1171       if(aNbExt > 0) {
1172         iLower=1;
1173         for (i=1; i<=aNbExt; i++) {
1174           aDist2=anExtrema.SquareDistance(i);
1175           if (aDist2 < aMinDist2) {
1176             aMinDist2=aDist2;
1177             iLower=i;
1178           }
1179         }
1180         aDist2=anExtrema.SquareDistance(iLower);
1181         Extrema_POnCurv aPOnC;
1182         Extrema_POnSurf aPOnS;
1183         anExtrema.Points(iLower, aPOnC, aPOnS);
1184         aTx=aPOnC.Parameter();
1185       }
1186       else {
1187         // modified by NIZHNY-MKK  Thu Jul 21 11:35:32 2005.BEGIN
1188         IntCurveSurface_HInter anExactIntersector;
1189
1190         Handle(GeomAdaptor_HCurve) aCurve     = new GeomAdaptor_HCurve(TheCurve);
1191         Handle(GeomAdaptor_HSurface) aSurface = new GeomAdaptor_HSurface(TheSurface);
1192
1193         anExactIntersector.Perform(aCurve, aSurface);
1194
1195         if(anExactIntersector.IsDone()) {
1196           Standard_Integer i = 0;
1197
1198           for(i = 1; i <= anExactIntersector.NbPoints(); i++) {
1199             const IntCurveSurface_IntersectionPoint& aPoint = anExactIntersector.Point(i);
1200
1201             if((aPoint.W() >= aTF) && (aPoint.W() <= aTL)) {
1202               aDist2=0.;
1203               aTx = aPoint.W();
1204             }
1205           }
1206         }
1207         // modified by NIZHNY-MKK  Thu Jul 21 11:35:40 2005.END
1208       }
1209     }
1210     else {
1211       return theflag;
1212     }
1213   }
1214
1215   Standard_Real aBoundaryDist;
1216
1217   aBoundaryDist = DistanceFunction(aTF) + myCriteria;
1218   if(aBoundaryDist * aBoundaryDist < aDist2) {
1219     aDist2 = aBoundaryDist * aBoundaryDist;
1220     aTx = aTF;
1221   }
1222
1223   aBoundaryDist = DistanceFunction(aTL) + myCriteria;
1224   if(aBoundaryDist * aBoundaryDist < aDist2) {
1225     aDist2 = aBoundaryDist * aBoundaryDist;
1226     aTx = aTL;
1227   }
1228
1229   Standard_Real aParameter = (aTF + aTL) * 0.5;
1230   aBoundaryDist = DistanceFunction(aParameter) + myCriteria;
1231   if(aBoundaryDist * aBoundaryDist < aDist2) {
1232     aDist2 = aBoundaryDist * aBoundaryDist;
1233     aTx = aParameter;
1234   }
1235
1236   if(aDist2 > myCriteria * myCriteria) {
1237     return theflag;
1238   }
1239
1240   if (fabs (aTx-aTF) < myEpsT) {
1241     return !theflag;
1242   }
1243
1244   if (fabs (aTx-aTL) < myEpsT) {
1245     return !theflag;
1246   }
1247
1248   if (aTx>aTF && aTx<aTL) {
1249     return !theflag;
1250   }
1251
1252   return theflag;
1253 }
1254
1255
1256 //=======================================================================
1257 //function : Perform
1258 //purpose  :
1259 //=======================================================================
1260   void IntTools_EdgeFace::Perform()
1261 {
1262   Standard_Integer i, aNb;
1263   IntTools_CommonPrt aCommonPrt;
1264   //
1265   aCommonPrt.SetEdge1(myEdge);
1266   //
1267   myErrorStatus=0;
1268   CheckData();
1269   if (myErrorStatus) {
1270     return;
1271   }
1272   //
1273   if (myContext.IsNull()) {
1274     myContext=new BOPInt_Context;
1275   }
1276   //
1277   myIsDone = Standard_False;
1278   myC.Initialize(myEdge);
1279   GeomAbs_CurveType aCurveType;
1280   aCurveType=myC.GetType();
1281   //
1282   // Prepare myCriteria
1283   if (aCurveType==GeomAbs_BSplineCurve||
1284         aCurveType==GeomAbs_BezierCurve) {
1285     //--- 5112
1286     Standard_Real diff1 = (myTolE/myTolF);
1287     Standard_Real diff2 = (myTolF/myTolE);
1288     if( diff1 > 100 || diff2 > 100 ) {
1289       myCriteria = Max(myTolE,myTolF);
1290     }
1291     else //--- 5112
1292       myCriteria=1.5*myTolE+myTolF;
1293   }
1294   else {
1295     myCriteria=myTolE+myTolF;
1296   }
1297
1298   myTmin=myRange.First();
1299   myTmax=myRange.Last();
1300
1301   myS.Initialize (myFace,Standard_True);
1302
1303   if(myContext.IsNull()) {
1304     myFClass2d.Init(myFace, 1.e-6);
1305   }
1306
1307   IntTools_BeanFaceIntersector anIntersector(myC, myS, myTolE, myTolF);
1308   anIntersector.SetBeanParameters(myRange.First(), myRange.Last());
1309   //
1310   anIntersector.SetContext(myContext);
1311   //
1312   anIntersector.Perform();
1313
1314   if(!anIntersector.IsDone()) {
1315     return;
1316   }
1317
1318   for(Standard_Integer r = 1; r <= anIntersector.Result().Length(); r++) {
1319     const IntTools_Range& aRange = anIntersector.Result().Value(r);
1320
1321     if(IsProjectable(IntTools_Tools::IntermediatePoint(aRange.First(), aRange.Last()))) {
1322       aCommonPrt.SetRange1(aRange.First(), aRange.Last());
1323       mySeqOfCommonPrts.Append(aCommonPrt);
1324     }
1325   }
1326
1327   aNb = mySeqOfCommonPrts.Length();
1328
1329   for (i=1; i<=aNb; i++) {
1330     IntTools_CommonPrt& aCP=mySeqOfCommonPrts.ChangeValue(i);
1331     //
1332     Standard_Real aTx1, aTx2;
1333     gp_Pnt aPx1, aPx2;
1334     //
1335     aCP.Range1(aTx1, aTx2);
1336     myC.D0(aTx1, aPx1);
1337     myC.D0(aTx2, aPx2);
1338     aCP.SetBoundingPoints(aPx1, aPx2);
1339     //
1340     MakeType (aCP);
1341   }
1342   {
1343     // Line\Cylinder's Common Parts treatement
1344     GeomAbs_CurveType   aCType;
1345     GeomAbs_SurfaceType aSType;
1346     TopAbs_ShapeEnum aType;
1347     Standard_Boolean bIsTouch;
1348     Standard_Real aTx;
1349
1350     aCType=myC.GetType();
1351     aSType=myS.GetType();
1352
1353     if (aCType==GeomAbs_Line && aSType==GeomAbs_Cylinder) {
1354       for (i=1; i<=aNb; i++) {
1355         IntTools_CommonPrt& aCP=mySeqOfCommonPrts(i);
1356         aType=aCP.Type();
1357         if (aType==TopAbs_EDGE) {
1358           bIsTouch=CheckTouch (aCP, aTx);
1359           if (bIsTouch) {
1360             aCP.SetType(TopAbs_VERTEX);
1361             aCP.SetVertexParameter1(aTx);
1362             aCP.SetRange1 (aTx, aTx);
1363           }
1364         }
1365         if (aType==TopAbs_VERTEX) {
1366           bIsTouch=CheckTouchVertex (aCP, aTx);
1367           if (bIsTouch) {
1368             aCP.SetVertexParameter1(aTx);
1369             aCP.SetRange1 (aTx, aTx);
1370           }
1371         }
1372       }
1373     }
1374
1375     // Circle\Plane's Common Parts treatement
1376
1377     if (aCType==GeomAbs_Circle && aSType==GeomAbs_Plane) {
1378       Standard_Boolean bIsCoplanar, bIsRadius;
1379       bIsCoplanar=IsCoplanar(myC, myS);
1380       bIsRadius=IsRadius(myC, myS);
1381       if (!bIsCoplanar && !bIsRadius) {
1382         for (i=1; i<=aNb; i++) {
1383           IntTools_CommonPrt& aCP=mySeqOfCommonPrts(i);
1384           aType=aCP.Type();
1385           if (aType==TopAbs_EDGE) {
1386             bIsTouch=CheckTouch (aCP, aTx);
1387             if (bIsTouch) {
1388               aCP.SetType(TopAbs_VERTEX);
1389               aCP.SetVertexParameter1(aTx);
1390               aCP.SetRange1 (aTx, aTx);
1391             }
1392           }
1393         }
1394       }
1395     }
1396   }
1397   myIsDone=Standard_True;
1398 }
1399
1400 //
1401 // myErrorStatus
1402 // 1 - the method Perform() is not invoked
1403 // 2,3,4,5 -the method CheckData() fails
1404 // 6 - PrepareArgs() problems
1405 // 7 - No Projectable ranges
1406 // 8,9 - PrepareArgs() problems occured inside  projectable Ranges
1407 // 11 - can't fill array  aFunc(i) in PrepareArgsFuncArrays
1408
1409
1410 //=======================================================================
1411 //function : CheckTouch
1412 //purpose  :
1413 //=======================================================================
1414   Standard_Boolean IntTools_EdgeFace::CheckTouchVertex (const IntTools_CommonPrt& aCP,
1415                                                         Standard_Real& aTx)
1416 {
1417   Standard_Real aTF, aTL, Tol, U1f,U1l,V1f,V1l, af, al,aDist2, aMinDist2, aTm, aDist2New;
1418   Standard_Real aEpsT;
1419   Standard_Boolean theflag=Standard_False;
1420   Standard_Integer aNbExt, i, iLower ;
1421
1422   aCP.Range1(aTF, aTL);
1423   aEpsT=8.e-5;
1424   aTm=0.5*(aTF+aTL);
1425   aDist2=DistanceFunction(aTm);
1426   aDist2 *= aDist2;
1427
1428   Tol = Precision::PConfusion();
1429
1430   const Handle(Geom_Curve)&  Curve   =BRep_Tool::Curve  (myC.Edge(), af, al);
1431   const Handle(Geom_Surface)& Surface=BRep_Tool::Surface(myS.Face());
1432
1433   Surface->Bounds(U1f,U1l,V1f,V1l);
1434
1435   GeomAdaptor_Curve   TheCurve   (Curve,aTF, aTL);
1436   GeomAdaptor_Surface TheSurface (Surface, U1f, U1l, V1f, V1l);
1437
1438   Extrema_ExtCS anExtrema (TheCurve, TheSurface, Tol, Tol);
1439
1440   if(!anExtrema.IsDone()) {
1441     return theflag;
1442   }
1443   if (anExtrema.IsParallel()) {
1444     return theflag;
1445   }
1446
1447   aNbExt=anExtrema.NbExt() ;
1448   if (!aNbExt) {
1449      return theflag;
1450   }
1451
1452   iLower=1;
1453   aMinDist2=1.e100;
1454   for (i=1; i<=aNbExt; ++i) {
1455     aDist2=anExtrema.SquareDistance(i);
1456     if (aDist2 < aMinDist2) {
1457       aMinDist2=aDist2;
1458       iLower=i;
1459     }
1460   }
1461
1462   aDist2New=anExtrema.SquareDistance(iLower);
1463
1464   if (aDist2New > aDist2) {
1465     aTx=aTm;
1466     return !theflag;
1467   }
1468
1469   if (aDist2New > myCriteria * myCriteria) {
1470     return theflag;
1471   }
1472
1473   Extrema_POnCurv aPOnC;
1474   Extrema_POnSurf aPOnS;
1475   anExtrema.Points(iLower, aPOnC, aPOnS);
1476
1477
1478   aTx=aPOnC.Parameter();
1479   ///
1480   if (fabs (aTx-aTF) < aEpsT) {
1481     return theflag;
1482   }
1483
1484   if (fabs (aTx-aTL) < aEpsT) {
1485     return theflag;
1486   }
1487
1488   if (aTx>aTF && aTx<aTL) {
1489     return !theflag;
1490   }
1491
1492   return theflag;
1493 }
1494
1495
1496 //=======================================================================
1497 //function :  IsCoplanar
1498 //purpose  :
1499 //=======================================================================
1500 Standard_Boolean IsCoplanar (const BRepAdaptor_Curve& aCurve ,
1501                                const BRepAdaptor_Surface& aSurface)
1502 {
1503   Standard_Boolean bFlag=Standard_False;
1504
1505   GeomAbs_CurveType   aCType;
1506   GeomAbs_SurfaceType aSType;
1507
1508   aCType=aCurve.GetType();
1509   aSType=aSurface.GetType();
1510
1511   if (aCType==GeomAbs_Circle && aSType==GeomAbs_Plane) {
1512     gp_Circ aCirc=aCurve.Circle();
1513     const gp_Ax1& anAx1=aCirc.Axis();
1514     const gp_Dir& aDirAx1=anAx1.Direction();
1515
1516     gp_Pln aPln=aSurface.Plane();
1517     const gp_Ax1& anAx=aPln.Axis();
1518     const gp_Dir& aDirPln=anAx.Direction();
1519
1520     bFlag=IntTools_Tools::IsDirsCoinside(aDirAx1, aDirPln);
1521   }
1522   return bFlag;
1523 }
1524 //=======================================================================
1525 //function :  IsRadius
1526 //purpose  :
1527 //=======================================================================
1528 Standard_Boolean IsRadius (const BRepAdaptor_Curve& aCurve ,
1529                            const BRepAdaptor_Surface& aSurface)
1530 {
1531   Standard_Boolean bFlag=Standard_False;
1532
1533   GeomAbs_CurveType   aCType;
1534   GeomAbs_SurfaceType aSType;
1535
1536   aCType=aCurve.GetType();
1537   aSType=aSurface.GetType();
1538
1539   if (aCType==GeomAbs_Circle && aSType==GeomAbs_Plane) {
1540     gp_Circ aCirc=aCurve.Circle();
1541     const gp_Pnt aCenter=aCirc.Location();
1542     Standard_Real aR=aCirc.Radius();
1543     gp_Pln aPln=aSurface.Plane();
1544     Standard_Real aD=aPln.Distance(aCenter);
1545     if (fabs (aD-aR) < 1.e-7) {
1546       return !bFlag;
1547     }
1548   }
1549   return bFlag;
1550 }
1551 //
1552 //=======================================================================
1553 //function :  AdaptiveDiscret
1554 //purpose  :
1555 //=======================================================================
1556 Standard_Integer AdaptiveDiscret (const Standard_Integer iDiscret,
1557                                   const BRepAdaptor_Curve& aCurve ,
1558                                   const BRepAdaptor_Surface& aSurface)
1559 {
1560   Standard_Integer iDiscretNew;
1561
1562   iDiscretNew=iDiscret;
1563
1564   GeomAbs_CurveType   aCType;
1565   GeomAbs_SurfaceType aSType;
1566
1567   aCType=aCurve.GetType();
1568   aSType=aSurface.GetType();
1569
1570   if (aSType==GeomAbs_Cylinder) {
1571    Standard_Real aELength, aRadius, dL, dLR;
1572
1573    aELength=IntTools::Length(aCurve.Edge());
1574    dL=aELength/iDiscret;
1575
1576    gp_Cylinder aCylinder=aSurface.Cylinder();
1577    aRadius=aCylinder.Radius();
1578    dLR=2*aRadius;
1579
1580    iDiscretNew=(Standard_Integer)(aELength/dLR);
1581
1582    if (iDiscretNew<iDiscret) {
1583      iDiscretNew=iDiscret;
1584    }
1585
1586   }
1587   return iDiscretNew;
1588 }
1589
1590
1591 #ifdef WNT
1592 #pragma warning ( default : 4101 )
1593 #endif