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