0026619: Tolerances of operands are modified using bop
[occt.git] / src / IntTools / IntTools_EdgeEdge.cxx
1 // Created by: Eugeny MALTCHIKOV
2 // Copyright (c) 2013-2014 OPEN CASCADE SAS
3 //
4 // This file is part of Open CASCADE Technology software library.
5 //
6 // This library is free software; you can redistribute it and/or modify it under
7 // the terms of the GNU Lesser General Public License version 2.1 as published
8 // by the Free Software Foundation, with special exception defined in the file
9 // OCCT_LGPL_EXCEPTION.txt. Consult the file LICENSE_LGPL_21.txt included in OCCT
10 // distribution for complete text of the license and disclaimer of any warranty.
11 //
12 // Alternatively, this file may be used under the terms of Open CASCADE
13 // commercial license or contractual agreement.
14
15
16 #include <Bnd_Box.hxx>
17 #include <BndLib_Add3dCurve.hxx>
18 #include <BOPCol_MapOfInteger.hxx>
19 #include <BRep_Tool.hxx>
20 #include <BRepAdaptor_Curve.hxx>
21 #include <ElCLib.hxx>
22 #include <Geom_BezierCurve.hxx>
23 #include <Geom_BSplineCurve.hxx>
24 #include <Geom_Circle.hxx>
25 #include <Geom_Curve.hxx>
26 #include <Geom_Ellipse.hxx>
27 #include <GeomAPI_ProjectPointOnCurve.hxx>
28 #include <gp_Dir.hxx>
29 #include <gp_Lin.hxx>
30 #include <IntTools_CommonPrt.hxx>
31 #include <IntTools_EdgeEdge.hxx>
32 #include <IntTools_Range.hxx>
33 #include <IntTools_Tools.hxx>
34 #include <TopoDS_Edge.hxx>
35 #include <TopoDS_Iterator.hxx>
36
37 static 
38   void BndBuildBox(const BRepAdaptor_Curve& theBAC,
39                    const Standard_Real aT1,
40                    const Standard_Real aT2,
41                    const Standard_Real theTol,
42                    Bnd_Box& theBox);
43 static
44   Standard_Real PointBoxDistance(const Bnd_Box& aB,
45                                  const gp_Pnt& aP);
46 static 
47   Standard_Integer SplitRangeOnSegments(const Standard_Real aT1, 
48                                         const Standard_Real aT2,
49                                         const Standard_Real theResolution,
50                                         const Standard_Integer theNbSeg,
51                                         IntTools_SequenceOfRanges& theSegments);
52 static
53  Standard_Integer DistPC(const Standard_Real aT1, 
54                          const Handle(Geom_Curve)& theC1,
55                          const Standard_Real theCriteria, 
56                          GeomAPI_ProjectPointOnCurve& theProjector,
57                          Standard_Real& aD, 
58                          Standard_Real& aT2,
59                          const Standard_Integer iC = 1);
60 static
61  Standard_Integer DistPC(const Standard_Real aT1, 
62                          const Handle(Geom_Curve)& theC1,
63                          const Standard_Real theCriteria,
64                          GeomAPI_ProjectPointOnCurve& theProjector, 
65                          Standard_Real& aD, 
66                          Standard_Real& aT2,
67                          Standard_Real& aDmax,
68                          Standard_Real& aT1max,
69                          Standard_Real& aT2max,
70                          const Standard_Integer iC = 1);
71 static
72   Standard_Integer FindDistPC(const Standard_Real aT1A, 
73                               const Standard_Real aT1B,
74                               const Handle(Geom_Curve)& theC1,
75                               const Standard_Real theCriteria,
76                               const Standard_Real theEps,
77                               GeomAPI_ProjectPointOnCurve& theProjector,
78                               Standard_Real& aDmax, 
79                               Standard_Real& aT1max,
80                               Standard_Real& aT2max,
81                               const Standard_Boolean bMaxDist = Standard_True);
82 static
83   Standard_Real ResolutionCoeff(const BRepAdaptor_Curve& theBAC,
84                                 const IntTools_Range& theRange);
85 static
86   Standard_Real Resolution(const Handle(Geom_Curve)& theCurve,
87                            const GeomAbs_CurveType theCurveType,
88                            const Standard_Real theResCoeff,
89                            const Standard_Real theR3D);
90 static
91   Standard_Real CurveDeflection(const BRepAdaptor_Curve& theBAC,
92                                 const IntTools_Range& theRange);
93 static 
94   Standard_Integer IsClosed(const Handle(Geom_Curve)& theCurve,
95                             const Standard_Real aT1,
96                             const Standard_Real aT2,
97                             const Standard_Real theTol,
98                             const Standard_Real theRes);
99 static 
100   Standard_Integer TypeToInteger(const GeomAbs_CurveType theCType);
101
102 //=======================================================================
103 //function : Prepare
104 //purpose  : 
105 //=======================================================================
106 void IntTools_EdgeEdge::Prepare()
107 {
108   GeomAbs_CurveType aCT1, aCT2;
109   Standard_Integer iCT1, iCT2;
110   //
111   myCurve1.Initialize(myEdge1);
112   myCurve2.Initialize(myEdge2);
113   //
114   if (myRange1.First() == 0. && myRange1.Last() == 0.) {
115     myRange1.SetFirst(myCurve1.FirstParameter());
116     myRange1.SetLast (myCurve1.LastParameter());
117   }
118   //
119   if (myRange2.First() == 0. && myRange2.Last() == 0.) {
120     myRange2.SetFirst(myCurve2.FirstParameter());
121     myRange2.SetLast (myCurve2.LastParameter());
122   }
123   //
124   aCT1 = myCurve1.GetType();
125   aCT2 = myCurve2.GetType();
126   //
127   iCT1 = TypeToInteger(aCT1);
128   iCT2 = TypeToInteger(aCT2);
129   //
130   if (iCT1 == iCT2) {
131     if (iCT1 != 0) {
132       //compute deflection
133       Standard_Real aC1, aC2;
134       //
135       aC2 = CurveDeflection(myCurve2, myRange2);
136       aC1 = (aC2 > Precision::Confusion()) ? 
137         CurveDeflection(myCurve1, myRange1) : 1.;
138       //
139       if (aC1 < aC2) {
140         --iCT1;
141       }
142     }
143   }
144   //
145   if (iCT1 < iCT2) {
146     TopoDS_Edge tmpE = myEdge1;
147     myEdge1 = myEdge2;
148     myEdge2 = tmpE;
149     //
150     BRepAdaptor_Curve tmpC = myCurve1;
151     myCurve1 = myCurve2;
152     myCurve2 = tmpC;
153     //
154     IntTools_Range tmpR = myRange1;
155     myRange1 = myRange2;
156     myRange2 = tmpR;
157     //
158     mySwap = Standard_True;
159   }
160   //
161   Standard_Real aTolAdd = Precision::Confusion() / 2.;
162   myTol1 = myCurve1.Tolerance() + aTolAdd;
163   myTol2 = myCurve2.Tolerance() + aTolAdd;
164   myTol = myTol1 + myTol2;
165   //
166   if (iCT1 != 0 || iCT2 != 0) {
167     Standard_Real f, l, aTM;
168     //
169     myGeom1 = BRep_Tool::Curve(myEdge1, f, l);
170     myGeom2 = BRep_Tool::Curve(myEdge2, f, l);
171     //
172     myResCoeff1 = ResolutionCoeff(myCurve1, myRange1);
173     myResCoeff2 = ResolutionCoeff(myCurve2, myRange2);
174     //
175     myRes1 = Resolution(myCurve1.Curve().Curve(), myCurve1.GetType(), myResCoeff1, myTol1);
176     myRes2 = Resolution(myCurve2.Curve().Curve(), myCurve2.GetType(), myResCoeff2, myTol2);
177     //
178     myPTol1 = 5.e-13;
179     aTM = Max(fabs(myRange1.First()), fabs(myRange1.Last()));
180     if (aTM > 999.) {
181       myPTol1 = 5.e-16 * aTM;
182     }
183     //
184     myPTol2 = 5.e-13;
185     aTM = Max(fabs(myRange2.First()), fabs(myRange2.Last()));
186     if (aTM > 999.) {
187       myPTol2 = 5.e-16 * aTM;
188     }
189   }
190 }
191
192 //=======================================================================
193 //function : Perform
194 //purpose  : 
195 //=======================================================================
196 void IntTools_EdgeEdge::Perform()
197 {
198   //1. Check data
199   CheckData();
200   if (myErrorStatus) {
201     return;
202   }
203   //
204   //2. Prepare Data
205   Prepare();
206   //
207   //3.1. Check Line/Line case
208   if (myCurve1.GetType() == GeomAbs_Line &&
209       myCurve2.GetType() == GeomAbs_Line) {
210     ComputeLineLine();
211     return;
212   }
213   //
214   IntTools_SequenceOfRanges aRanges1, aRanges2;
215   //
216   //3.2. Find ranges containig solutions
217   Standard_Boolean bSplit2;
218   FindSolutions(aRanges1, aRanges2, bSplit2);
219   //
220   //4. Merge solutions and save common parts
221   MergeSolutions(aRanges1, aRanges2, bSplit2);
222 }
223
224 //=======================================================================
225 //function : FindSolutions
226 //purpose  : 
227 //=======================================================================
228 void IntTools_EdgeEdge::FindSolutions(IntTools_SequenceOfRanges& theRanges1,
229                                       IntTools_SequenceOfRanges& theRanges2,
230                                       Standard_Boolean& bSplit2)
231 {
232   Standard_Boolean bIsClosed2;
233   Standard_Real aT11, aT12, aT21, aT22;
234   Bnd_Box aB2;
235   //
236   bSplit2 = Standard_False;
237   myRange1.Range(aT11, aT12);
238   myRange2.Range(aT21, aT22);
239   //
240   bIsClosed2 = IsClosed(myGeom2, aT21, aT22, myTol2, myRes2);
241   //
242   if (bIsClosed2) {
243     Bnd_Box aB1;
244     BndBuildBox(myCurve1, aT11, aT12, myTol1, aB1);
245     //
246     gp_Pnt aP = myGeom2->Value(aT21);
247     bIsClosed2 = !aB1.IsOut(aP);
248   }
249   //
250   if (!bIsClosed2) {
251     BndBuildBox(myCurve2, aT21, aT22, myTol2, aB2);
252     FindSolutions(myRange1, myRange2, aB2, theRanges1, theRanges2);
253     return;
254   }
255   //
256   if (!CheckCoincidence(aT11, aT12, aT21, aT22, myTol, myRes1)) {
257     theRanges1.Append(myRange1);
258     theRanges2.Append(myRange2);
259     return;
260   }
261   //
262   Standard_Integer i, j, aNb1, aNb2;
263   IntTools_SequenceOfRanges aSegments1, aSegments2;
264   //
265   aNb1 = IsClosed(myGeom1, aT11, aT12, myTol1, myRes1) ? 2 : 1;
266   aNb2 = 2;
267   //
268   aNb1 = SplitRangeOnSegments(aT11, aT12, myRes1, aNb1, aSegments1);
269   aNb2 = SplitRangeOnSegments(aT21, aT22, myRes2, aNb2, aSegments2);
270   //
271   for (i = 1; i <= aNb1; ++i) {
272     const IntTools_Range& aR1 = aSegments1(i);
273     for (j = 1; j <= aNb2; ++j) {
274       const IntTools_Range& aR2 = aSegments2(j);
275       BndBuildBox(myCurve2, aR2.First(), aR2.Last(), myTol2, aB2);
276       FindSolutions(aR1, aR2, aB2, theRanges1, theRanges2);
277     }
278   }
279   //
280   bSplit2 = aNb2 > 1;
281 }
282
283 //=======================================================================
284 //function : FindSolutions
285 //purpose  : 
286 //=======================================================================
287 void IntTools_EdgeEdge::FindSolutions(const IntTools_Range& theR1,
288                                       const IntTools_Range& theR2,
289                                       const Bnd_Box& theBox2,
290                                       IntTools_SequenceOfRanges& theRanges1,
291                                       IntTools_SequenceOfRanges& theRanges2)
292 {
293   Standard_Boolean bOut, bStop, bThin;
294   Standard_Real aT11, aT12, aT21, aT22;
295   Standard_Real aTB11, aTB12, aTB21, aTB22;
296   Standard_Real aSmallStep1, aSmallStep2;
297   Standard_Integer iCom;
298   Bnd_Box aB1, aB2;
299   //
300   theR1.Range(aT11, aT12);
301   theR2.Range(aT21, aT22);
302   //
303   aB2 = theBox2;
304   //
305   bThin = Standard_False;
306   bStop = Standard_False;
307   iCom  = 1;
308   //
309   do {
310     aTB11 = aT11;
311     aTB12 = aT12;
312     aTB21 = aT21;
313     aTB22 = aT22;
314     //
315     //1. Build box for first edge and find parameters 
316     //   of the second one in that box
317     BndBuildBox(myCurve1, aT11, aT12, myTol1, aB1);
318     bOut = aB1.IsOut(aB2);
319     if (bOut) {
320       break;
321     }
322     //
323     bThin = ((aT12 - aT11) < myRes1) ||
324       (aB1.IsXThin(myTol) && aB1.IsYThin(myTol) && aB1.IsZThin(myTol));
325     //
326     bOut = !FindParameters(myCurve2, aTB21, aTB22, myTol2, myRes2, myPTol2, 
327                            myResCoeff2, aB1, aT21, aT22);
328     if (bOut || bThin) {
329       break;
330     }
331     //
332     //2. Build box for second edge and find parameters 
333     //   of the first one in that box
334     BndBuildBox(myCurve2, aT21, aT22, myTol2, aB2);
335     bOut = aB1.IsOut(aB2);
336     if (bOut) {
337       break;
338     }
339     //
340     bThin = ((aT22 - aT21) < myRes2) ||
341       (aB2.IsXThin(myTol) && aB2.IsYThin(myTol) && aB2.IsZThin(myTol));
342     //
343     bOut = !FindParameters(myCurve1, aTB11, aTB12, myTol1, myRes1, myPTol1,
344                            myResCoeff1, aB2, aT11, aT12);
345     //
346     if (bOut || bThin) {
347       break;
348     }
349     //
350     //3. Check if it makes sense to continue
351     aSmallStep1 = (aTB12 - aTB11) / 250.;
352     aSmallStep2 = (aTB22 - aTB21) / 250.;
353     //
354     if (aSmallStep1 < myRes1) {
355       aSmallStep1 = myRes1;
356     }
357     if (aSmallStep2 < myRes2) {
358       aSmallStep2 = myRes2;
359     }
360     //
361     if (((aT11 - aTB11) < aSmallStep1) && ((aTB12 - aT12) < aSmallStep1) &&
362         ((aT21 - aTB21) < aSmallStep2) && ((aTB22 - aT22) < aSmallStep2)) {
363       bStop = Standard_True;
364     }
365     //
366   } while (!bStop);
367   //
368   if (bOut) {
369     //no intersection;
370     return;
371   }
372   //
373   if (!bThin) {
374     //check curves for coincidence on the ranges
375     iCom = CheckCoincidence(aT11, aT12, aT21, aT22, myTol, myRes1);
376     if (!iCom) {
377       bThin = Standard_True;
378     }
379   }
380   //
381   if (bThin) {
382     if (iCom != 0) {
383       //check intermediate points
384       Standard_Boolean bSol;
385       Standard_Real aT1;
386       gp_Pnt aP1;
387       GeomAPI_ProjectPointOnCurve aProjPC;
388       //
389       aT1 = (aT11 + aT12) * .5;
390       myGeom1->D0(aT1, aP1);
391       //
392       aProjPC.Init(myGeom2, aT21, aT22);
393       aProjPC.Perform(aP1);
394       //
395       if (aProjPC.NbPoints()) {
396         bSol = aProjPC.LowerDistance() <= myTol;
397       }
398       else {
399         Standard_Real aT2;
400         gp_Pnt aP2;
401         //
402         aT2 = (aT21 + aT22) * .5;
403         myGeom2->D0(aT2, aP2);
404         //
405         bSol = aP1.IsEqual(aP2, myTol);
406       }
407       //
408       if (!bSol) {
409         return;
410       }
411     }
412     //add common part
413     IntTools_Range aR1(aT11, aT12), aR2(aT21, aT22);
414     //
415     theRanges1.Append(aR1);
416     theRanges2.Append(aR2);
417     return;
418   }
419   //
420   if (!IsIntersection(aT11, aT12, aT21, aT22)) {
421     return;
422   }
423   //
424   //split ranges on segments and repeat
425   Standard_Integer i, aNb1;
426   IntTools_SequenceOfRanges aSegments1;
427   //
428   IntTools_Range aR2(aT21, aT22);
429   BndBuildBox(myCurve2, aT21, aT22, myTol2, aB2);
430   //
431   aNb1 = SplitRangeOnSegments(aT11, aT12, myRes1, 3, aSegments1);
432   for (i = 1; i <= aNb1; ++i) {
433     const IntTools_Range& aR1 = aSegments1(i);
434     FindSolutions(aR1, aR2, aB2, theRanges1, theRanges2);
435   }
436 }
437
438 //=======================================================================
439 //function : FindParameters
440 //purpose  : 
441 //=======================================================================
442 Standard_Boolean IntTools_EdgeEdge::FindParameters(const BRepAdaptor_Curve& theBAC,
443                                                    const Standard_Real aT1, 
444                                                    const Standard_Real aT2,
445                                                    const Standard_Real theTol,
446                                                    const Standard_Real theRes,
447                                                    const Standard_Real thePTol,
448                                                    const Standard_Real theResCoeff,
449                                                    const Bnd_Box& theCBox,
450                                                    Standard_Real& aTB1, 
451                                                    Standard_Real& aTB2)
452 {
453   Standard_Boolean bRet;
454   Standard_Integer aC, i, k;
455   Standard_Real aCf, aDiff, aDt, aT, aTB, aTOut, aTIn;
456   Standard_Real aDist, aDistP, aDistTol;
457   gp_Pnt aP;
458   Bnd_Box aCBx;
459   //
460   bRet = Standard_False;
461   aCf = 0.6180339887498948482045868343656;// =0.5*(1.+sqrt(5.))/2.;
462   aDt = theRes;
463   aDistP = 0.;
464   aDistTol = 1e-9;
465   aCBx = theCBox;
466   aCBx.Enlarge(theTol);
467   //
468   const Handle(Geom_Curve)& aCurve = theBAC.Curve().Curve();
469   const GeomAbs_CurveType aCurveType = theBAC.GetType();
470   //
471   for (i = 0; i < 2; ++i) {
472     aTB = !i ? aT1 : aT2;
473     aT = !i ? aT2 : aTB1;
474     aC = !i ? 1 : -1;
475     bRet = Standard_False;
476     k = 0;
477     //looking for the point on the edge which is in the box;
478     while (aC*(aT-aTB) >= 0) {
479       theBAC.D0(aTB, aP);
480       aDist = PointBoxDistance(theCBox, aP);
481       if (aDist > theTol) {
482         if (fabs(aDist - aDistP) < aDistTol) {
483           aDt = Resolution(aCurve, aCurveType, theResCoeff, (++k)*aDist);
484         } else {
485           k = 0;
486           aDt = Resolution(aCurve, aCurveType, theResCoeff, aDist);
487         }
488         aTB += aC*aDt;
489       } else {
490         bRet = Standard_True;
491         break;
492       }
493       aDistP = aDist;
494     }
495     //
496     if (!bRet) {
497       if (!i) {
498         //edge is out of the box;
499         return bRet;
500       } else {
501         bRet = !bRet;
502         aTB = aTB1;
503         aDt = aT2 - aTB1;
504       }
505     }
506     //
507     aT = !i ? aT1 : aT2;
508     if (aTB != aT) {
509       //one point IN, one point OUT; looking for the bounding point;
510       aTIn = aTB;
511       aTOut = aTB - aC*aDt;
512       aDiff = aTIn - aTOut;
513       while (fabs(aDiff) > thePTol) {
514         aTB = aTOut + aDiff*aCf;
515         theBAC.D0(aTB, aP);
516         if (aCBx.IsOut(aP)) {
517           aTOut = aTB;
518         } else {
519           aTIn = aTB;
520         }
521         aDiff = aTIn - aTOut;
522       }
523     }
524     if (!i) {
525       aTB1 = aTB;
526     } else {
527       aTB2 = aTB;
528     }
529   }
530   return bRet;
531 }
532
533 //=======================================================================
534 //function : MergeSolutions
535 //purpose  : 
536 //=======================================================================
537 void IntTools_EdgeEdge::MergeSolutions(const IntTools_SequenceOfRanges& theRanges1, 
538                                        const IntTools_SequenceOfRanges& theRanges2,
539                                        const Standard_Boolean bSplit2)
540 {
541   Standard_Integer aNbCP = theRanges1.Length();
542   if (aNbCP == 0) {
543     return;
544   }
545   //
546   IntTools_Range aRi1, aRi2, aRj1, aRj2;
547   Standard_Boolean bCond;
548   Standard_Integer i, j;
549   TopAbs_ShapeEnum aType;
550   Standard_Real aT11, aT12, aT21, aT22;
551   Standard_Real aTi11, aTi12, aTi21, aTi22;
552   Standard_Real aTj11, aTj12, aTj21, aTj22;
553   Standard_Real aRes1, aRes2, dTR1, dTR2;
554   BOPCol_MapOfInteger aMI;
555   //
556   aRes1 = Resolution(myCurve1.Curve().Curve(), 
557                      myCurve1.GetType(), myResCoeff1, myTol);
558   aRes2 = Resolution(myCurve2.Curve().Curve(), 
559                      myCurve2.GetType(), myResCoeff2, myTol);
560   //
561   myRange1.Range(aT11, aT12);
562   myRange2.Range(aT21, aT22);
563   dTR1 = 20*aRes1;
564   dTR2 = 20*aRes2;
565   aType = TopAbs_VERTEX;
566   //
567   for (i = 1; i <= aNbCP;) {
568     if (aMI.Contains(i)) {
569       ++i;
570       continue;
571     }
572     //
573     aRi1 = theRanges1(i);
574     aRi2 = theRanges2(i);
575     //
576     aRi1.Range(aTi11, aTi12);
577     aRi2.Range(aTi21, aTi22);
578     //
579     aMI.Add(i);
580     //
581     for (j = i+1; j <= aNbCP; ++j) {
582       if (aMI.Contains(j)) {
583         continue;
584       }
585       //
586       aRj1 = theRanges1(j);
587       aRj2 = theRanges2(j);
588       //
589       aRj1.Range(aTj11, aTj12);
590       aRj2.Range(aTj21, aTj22);
591       //
592       bCond = (fabs(aTi12 - aTj11) < dTR1) ||
593         (bSplit2 && (fabs(aTj12 - aTi11) < dTR1));
594       if (bCond && bSplit2) {
595         bCond = (fabs((Max(aTi22, aTj22) - Min(aTi21, aTj21)) - 
596                       ((aTi22 - aTi21) + (aTj22 - aTj21))) < dTR2);
597       }
598       //
599       if (bCond) {
600         aTi11 = Min(aTi11, aTj11);
601         aTi12 = Max(aTi12, aTj12);
602         aTi21 = Min(aTi21, aTj21);
603         aTi22 = Max(aTi22, aTj22);
604         aMI.Add(j);
605       }
606       else if (!bSplit2) {
607         i = j;
608         break;
609       }
610     }
611     //
612     if (((fabs(aT11 - aTi11) < myRes1) && (fabs(aT12 - aTi12) < myRes1)) ||
613         ((fabs(aT21 - aTi21) < myRes2) && (fabs(aT22 - aTi22) < myRes2))) {
614       aType = TopAbs_EDGE;
615       myCommonParts.Clear();
616     }
617     //
618     AddSolution(aTi11, aTi12, aTi21, aTi22, aType);
619     if (aType == TopAbs_EDGE) {
620       break;
621     }
622     //
623     if (bSplit2) {
624       ++i;
625     }
626   }
627 }
628
629 //=======================================================================
630 //function : AddSolution
631 //purpose  : 
632 //=======================================================================
633 void IntTools_EdgeEdge::AddSolution(const Standard_Real aT11,
634                                     const Standard_Real aT12,
635                                     const Standard_Real aT21,
636                                     const Standard_Real aT22,
637                                     const TopAbs_ShapeEnum theType)
638 {
639   IntTools_CommonPrt aCPart;
640   //
641   aCPart.SetType(theType);
642   if (!mySwap) {
643     aCPart.SetEdge1(myEdge1);
644     aCPart.SetEdge2(myEdge2);
645     aCPart.SetRange1(aT11, aT12);
646     aCPart.AppendRange2(aT21, aT22);
647   } else {
648     aCPart.SetEdge1(myEdge2);
649     aCPart.SetEdge2(myEdge1);
650     aCPart.SetRange1(aT21, aT22);
651     aCPart.AppendRange2(aT11, aT12);
652   }
653   //
654   if (theType == TopAbs_VERTEX) {
655     Standard_Real aT1, aT2;
656     //
657     FindBestSolution(aT11, aT12, aT21, aT22, aT1, aT2);
658     //
659     if (!mySwap) {
660       aCPart.SetVertexParameter1(aT1);
661       aCPart.SetVertexParameter2(aT2);
662     } else {
663       aCPart.SetVertexParameter1(aT2);
664       aCPart.SetVertexParameter2(aT1);
665     }
666   }
667   myCommonParts.Append(aCPart);
668 }
669
670 //=======================================================================
671 //function : FindBestSolution
672 //purpose  : 
673 //=======================================================================
674 void IntTools_EdgeEdge::FindBestSolution(const Standard_Real aT11,
675                                          const Standard_Real aT12,
676                                          const Standard_Real aT21,
677                                          const Standard_Real aT22,
678                                          Standard_Real& aT1,
679                                          Standard_Real& aT2)
680 {
681   Standard_Integer i, aNbS, iErr;
682   Standard_Real aDMin, aD, aRes1, aSolCriteria, aTouchCriteria;
683   Standard_Real aT1A, aT1B, aT1Min, aT2Min;
684   Standard_Real aT1Im, aT2Im, aT1Touch;
685   GeomAPI_ProjectPointOnCurve aProjPC;
686   IntTools_SequenceOfRanges aRanges;
687   Standard_Boolean bTouch;
688   //
689   aDMin = Precision::Infinite();
690   aSolCriteria   = 5.e-16;
691   aTouchCriteria = 5.e-13;
692   bTouch = Standard_False;
693   aT1Touch = aT11;
694   //
695   aRes1 = Resolution(myCurve1.Curve().Curve(), 
696                      myCurve1.GetType(), myResCoeff1, myTol);
697   aNbS = 10;
698   aNbS = SplitRangeOnSegments(aT11, aT12, 3*aRes1, aNbS, aRanges);
699   //
700   aProjPC.Init(myGeom2, aT21, aT22);
701   //
702   aT1 = (aT11 + aT12) * 0.5;
703   iErr = DistPC(aT1, myGeom1, aSolCriteria, aProjPC, aD, aT2, -1);
704   if (iErr == 1) {
705     aT2 = (aT21 + aT22) * 0.5;
706   }
707   //
708   aT1Im = aT1;
709   aT2Im = aT2;
710   //
711   for (i = 1; i <= aNbS; ++i) {
712     const IntTools_Range& aR1 = aRanges(i);
713     aR1.Range(aT1A, aT1B);
714     //
715     aD = myTol;
716     iErr = FindDistPC(aT1A, aT1B, myGeom1, aSolCriteria, myPTol1,
717                       aProjPC, aD, aT1Min, aT2Min, Standard_False);
718     if (iErr != 1) {
719       if (aD < aDMin) {
720         aT1 = aT1Min;
721         aT2 = aT2Min;
722         aDMin = aD;
723       }
724       //
725       if (aD < aTouchCriteria) {
726         if (bTouch) {
727           aT1A = (aT1Touch + aT1Min) * 0.5;
728           iErr = DistPC(aT1A, myGeom1, aTouchCriteria, 
729                         aProjPC, aD, aT2Min, -1);
730           if (aD > aTouchCriteria) {
731             aT1 = aT1Im;
732             aT2 = aT2Im;
733             break;
734           }
735         }
736         else {
737           aT1Touch = aT1Min;
738           bTouch = Standard_True;
739         }
740       }
741     }
742   }
743 }
744
745 //=======================================================================
746 //function : ComputeLineLine
747 //purpose  : 
748 //=======================================================================
749 void IntTools_EdgeEdge::ComputeLineLine()
750 {
751   Standard_Boolean IsParallel, IsCoincide;
752   Standard_Real aSin, aCos, aAng, aTol;
753   Standard_Real aT1, aT2, aT11, aT12, aT21, aT22;
754   gp_Pnt aP11, aP12;
755   gp_Lin aL1, aL2;
756   gp_Dir aD1, aD2;
757   IntTools_CommonPrt aCommonPrt;
758   //
759   IsParallel = Standard_False;
760   IsCoincide = Standard_False;
761   aTol = myTol*myTol;
762   aL1 = myCurve1.Line();
763   aL2 = myCurve2.Line();
764   aD1 = aL1.Position().Direction();
765   aD2 = aL2.Position().Direction();
766   myRange1.Range(aT11, aT12);
767   myRange2.Range(aT21, aT22);
768   //
769   aCommonPrt.SetEdge1(myEdge1);
770   aCommonPrt.SetEdge2(myEdge2);
771   //
772   aCos = aD1.Dot(aD2);
773   aAng = (aCos >= 0.) ? 2.*(1. - aCos) : 2.*(1. + aCos);
774   //
775   if(aAng <= Precision::Angular()) {
776     IsParallel = Standard_True;
777     if(aL1.SquareDistance(aL2.Location()) <= aTol) {
778       IsCoincide = Standard_True;
779       aP11 = ElCLib::Value(aT11, aL1);
780       aP12 = ElCLib::Value(aT12, aL1);
781     }
782   }
783   else {
784     aP11 = ElCLib::Value(aT11, aL1);
785     aP12 = ElCLib::Value(aT12, aL1);
786     if(aL2.SquareDistance(aP11) <= aTol && aL2.SquareDistance(aP12) <= aTol) {
787       IsCoincide = Standard_True;
788     }
789   }
790   //
791   if (IsCoincide) {
792     Standard_Real t21, t22;
793     //
794     t21 = ElCLib::Parameter(aL2, aP11);
795     t22 = ElCLib::Parameter(aL2, aP12);
796     if((t21 > aT22 && t22 > aT22) || (t21 < aT21 && t22 < aT21)) {
797       return;
798     }
799     //
800     Standard_Real temp;
801     if(t21 > t22) {
802       temp = t21;
803       t21 = t22;
804       t22 = temp;
805     }
806     //
807     if(t21 >= aT21) {
808       if(t22 <= aT22) {
809         aCommonPrt.SetRange1(aT11, aT12);
810         aCommonPrt.SetAllNullFlag(Standard_True);
811         aCommonPrt.AppendRange2(t21, t22);
812       }
813       else {
814         aCommonPrt.SetRange1(aT11, aT12 - (t22 - aT22));
815         aCommonPrt.AppendRange2(t21, aT22);
816       }
817     }
818     else {
819       aCommonPrt.SetRange1(aT11 + (aT21 - t21), aT12);
820       aCommonPrt.AppendRange2(aT21, t22);
821     }
822     aCommonPrt.SetType(TopAbs_EDGE);  
823     myCommonParts.Append(aCommonPrt);
824     return;
825   }
826   //
827   if (IsParallel) {
828     return;
829   }
830   //
831   {
832     TopoDS_Iterator aIt1, aIt2;
833     aIt1.Initialize(myEdge1);
834     for (; aIt1.More(); aIt1.Next()) {
835       const TopoDS_Shape& aV1 = aIt1.Value();
836       aIt2.Initialize(myEdge2);
837       for (; aIt2.More(); aIt2.Next()) {
838         const TopoDS_Shape& aV2 = aIt2.Value();
839         if (aV2.IsSame(aV1)) {
840           return;
841         }
842       }
843     }
844   }
845   //
846   aSin = 1. - aCos*aCos;
847   gp_Pnt O1 = aL1.Location();
848   gp_Pnt O2 = aL2.Location();
849   gp_Vec O1O2 (O1, O2);
850   //
851   aT2 = (aD1.XYZ()*(O1O2.Dot(aD1))-(O1O2.XYZ())).Dot(aD2.XYZ());
852   aT2 /= aSin;
853   //
854   if(aT2 < aT21 || aT2 > aT22) {
855     return;
856   }
857   //
858   gp_Pnt aP2(ElCLib::Value(aT2, aL2));
859   aT1 = (gp_Vec(O1, aP2)).Dot(aD1);
860   //
861   if(aT1 < aT11 || aT1 > aT12) {
862     return;
863   }
864   //
865   gp_Pnt aP1(ElCLib::Value(aT1, aL1));
866   Standard_Real aDist = aP1.SquareDistance(aP2);
867   //
868   if (aDist > aTol) {
869     return;
870   }
871   //
872   // compute correct range on the edges
873   Standard_Real anAngle, aDt1, aDt2;
874   //
875   anAngle = aD1.Angle(aD2);
876   //
877   aDt1 = IntTools_Tools::ComputeIntRange(myTol1, myTol2, anAngle);
878   aDt2 = IntTools_Tools::ComputeIntRange(myTol2, myTol1, anAngle);
879   //
880   aCommonPrt.SetRange1(aT1 - aDt1, aT1 + aDt1);
881   aCommonPrt.AppendRange2(aT2 - aDt2, aT2 + aDt2);
882   aCommonPrt.SetType(TopAbs_VERTEX);
883   aCommonPrt.SetVertexParameter1(aT1);
884   aCommonPrt.SetVertexParameter2(aT2);
885   myCommonParts.Append(aCommonPrt);
886 }
887
888 //=======================================================================
889 //function : IsIntersection
890 //purpose  : 
891 //=======================================================================
892 Standard_Boolean IntTools_EdgeEdge::IsIntersection(const Standard_Real aT11,
893                                                    const Standard_Real aT12,
894                                                    const Standard_Real aT21,
895                                                    const Standard_Real aT22)
896 {
897   Standard_Boolean bRet;
898   gp_Pnt aP11, aP12, aP21, aP22;
899   gp_Vec aV11, aV12, aV21, aV22;
900   Standard_Real aD11_21, aD11_22, aD12_21, aD12_22, aCriteria, aCoef;
901   Standard_Boolean bSmall_11_21, bSmall_11_22, bSmall_12_21, bSmall_12_22;
902   //
903   bRet = Standard_True;
904   aCoef = 1.e+5;
905   if (((aT12 - aT11) > aCoef*myRes1) && ((aT22 - aT21) > aCoef*myRes2)) {
906     aCoef = 5000;
907   } else {
908     Standard_Real aTRMin = Min((aT12 - aT11)/myRes1, (aT22 - aT21)/myRes2);
909     aCoef = aTRMin / 100.;
910     if (aCoef < 1.) {
911       aCoef = 1.;
912     }
913   }
914   aCriteria = aCoef * myTol;
915   aCriteria *= aCriteria;
916   //
917   myGeom1->D1(aT11, aP11, aV11);
918   myGeom1->D1(aT12, aP12, aV12);
919   myGeom2->D1(aT21, aP21, aV21);
920   myGeom2->D1(aT22, aP22, aV22);
921   //
922   aD11_21 = aP11.SquareDistance(aP21);
923   aD11_22 = aP11.SquareDistance(aP22);
924   aD12_21 = aP12.SquareDistance(aP21);
925   aD12_22 = aP12.SquareDistance(aP22);
926   //
927   bSmall_11_21 = aD11_21 < aCriteria;
928   bSmall_11_22 = aD11_22 < aCriteria;
929   bSmall_12_21 = aD12_21 < aCriteria;
930   bSmall_12_22 = aD12_22 < aCriteria;
931   //
932   if ((bSmall_11_21 && bSmall_12_22) ||
933       (bSmall_11_22 && bSmall_12_21)) {
934     if (aCoef == 1.) {
935       return bRet;
936     }
937     //
938     Standard_Real anAngleCriteria;
939     Standard_Real anAngle1, anAngle2;
940     //
941     anAngleCriteria = 5.e-3;
942     if (bSmall_11_21 && bSmall_12_22) {
943       anAngle1 = aV11.Angle(aV21);
944       anAngle2 = aV12.Angle(aV22);
945     } else {
946       anAngle1 = aV11.Angle(aV22);
947       anAngle2 = aV12.Angle(aV21);
948     }
949     //
950     if (((anAngle1 < anAngleCriteria) || ((M_PI - anAngle1) < anAngleCriteria)) ||
951         ((anAngle2 < anAngleCriteria) || ((M_PI - anAngle2) < anAngleCriteria))) {
952       GeomAPI_ProjectPointOnCurve aProjPC;
953       Standard_Integer iErr;
954       Standard_Real aD, aT1Min, aT2Min;
955       //
956       aD = Precision::Infinite();
957       aProjPC.Init(myGeom2, aT21, aT22);
958       iErr = FindDistPC(aT11, aT12, myGeom1, myTol, myRes1, 
959                         aProjPC, aD, aT1Min, aT2Min, Standard_False);
960       bRet = (iErr == 2);
961     }
962   }
963   return bRet;
964 }
965
966 //=======================================================================
967 //function : CheckCoincidence
968 //purpose  : 
969 //=======================================================================
970 Standard_Integer IntTools_EdgeEdge::CheckCoincidence(const Standard_Real aT11,
971                                                      const Standard_Real aT12,
972                                                      const Standard_Real aT21,
973                                                      const Standard_Real aT22,
974                                                      const Standard_Real theCriteria,
975                                                      const Standard_Real theCurveRes1)
976 {
977   Standard_Integer iErr, aNb, aNb1, i;
978   Standard_Real aT1A, aT1B, aT1max, aT2max, aDmax;
979   GeomAPI_ProjectPointOnCurve aProjPC;
980   IntTools_SequenceOfRanges aRanges;
981   //
982   iErr  = 0;
983   aDmax = -1.;
984   aProjPC.Init(myGeom2, aT21, aT22);
985   //
986   // 1. Express evaluation
987   aNb = 10; // Number of intervals on the curve #1
988   aNb1 = SplitRangeOnSegments(aT11, aT12, theCurveRes1, aNb, aRanges);
989   for (i = 1; i < aNb1; ++i) {
990     const IntTools_Range& aR1 = aRanges(i);
991     aR1.Range(aT1A, aT1B);
992     //
993     iErr = DistPC(aT1B, myGeom1, theCriteria, aProjPC, aDmax, aT2max);
994     if (iErr) {
995       return iErr;
996     }
997   }
998   //
999   // if the ranges in aRanges are less than theCurveRes1,
1000   // there is no need to do step 2 (deep evaluation)
1001   if (aNb1 < aNb) {
1002     return iErr;
1003   }
1004   //
1005   // 2. Deep evaluation
1006   for (i = 2; i < aNb1; ++i) {
1007     const IntTools_Range& aR1 = aRanges(i);
1008     aR1.Range(aT1A, aT1B);
1009     //
1010     iErr = FindDistPC(aT1A, aT1B, myGeom1, theCriteria, theCurveRes1, 
1011                       aProjPC, aDmax, aT1max, aT2max);
1012     if (iErr) {
1013       return iErr;
1014     }
1015   }
1016   // Possible values:
1017   // iErr == 0 - the patches are coincided
1018   // iErr == 1 - a point from aC1 can not be projected on aC2
1019   // iErr == 2 - the distance is too big
1020   return iErr;
1021 }
1022
1023 //=======================================================================
1024 //function : FindDistPC
1025 //purpose  : 
1026 //=======================================================================
1027 Standard_Integer FindDistPC(const Standard_Real aT1A, 
1028                             const Standard_Real aT1B,
1029                             const Handle(Geom_Curve)& theC1,
1030                             const Standard_Real theCriteria,
1031                             const Standard_Real theEps,
1032                             GeomAPI_ProjectPointOnCurve& theProjPC,
1033                             Standard_Real& aDmax, 
1034                             Standard_Real& aT1max,
1035                             Standard_Real& aT2max,
1036                             const Standard_Boolean bMaxDist) 
1037 {
1038   Standard_Integer iErr, iC;
1039   Standard_Real aGS, aXP, aA, aB, aXL, aYP, aYL, aT2P, aT2L;
1040   //
1041   iC = bMaxDist ? 1 : -1;
1042   iErr = 0;
1043   aT1max = aT2max = 0.; // silence GCC warning
1044   //
1045   aGS = 0.6180339887498948482045868343656;// =0.5*(1.+sqrt(5.))-1.;
1046   aA = aT1A;
1047   aB = aT1B;
1048   //
1049   // check bounds
1050   iErr = DistPC(aA, theC1, theCriteria, theProjPC, 
1051                 aYP, aT2P, aDmax, aT1max, aT2max, iC);
1052   if (iErr == 2) {
1053     return iErr;
1054   }
1055   //
1056   iErr = DistPC(aB, theC1, theCriteria, theProjPC, 
1057                 aYL, aT2L, aDmax, aT1max, aT2max, iC);
1058   if (iErr == 2) {
1059     return iErr;
1060   }
1061   //
1062   aXP = aA + (aB - aA)*aGS;
1063   aXL = aB - (aB - aA)*aGS;
1064   //
1065   iErr = DistPC(aXP, theC1, theCriteria, theProjPC, 
1066                 aYP, aT2P, aDmax, aT1max, aT2max, iC);
1067   if (iErr) {
1068     return iErr;
1069   }
1070   //
1071   iErr = DistPC(aXL, theC1, theCriteria, theProjPC, 
1072                 aYL, aT2L, aDmax, aT1max, aT2max, iC);
1073   if (iErr) {
1074     return iErr;
1075   }
1076   //
1077   for (;;) {
1078     if (iC*(aYP - aYL) > 0) {
1079       aA = aXL;
1080       aXL = aXP;
1081       aYL = aYP;
1082       aXP = aA + (aB - aA)*aGS;
1083       iErr = DistPC(aXP, theC1, theCriteria, theProjPC, 
1084                     aYP, aT2P, aDmax, aT1max, aT2max, iC);
1085     }
1086     else {
1087       aB = aXP;
1088       aXP = aXL;
1089       aYP = aYL;
1090       aXL = aB - (aB - aA)*aGS;
1091       iErr = DistPC(aXL, theC1, theCriteria, theProjPC, 
1092                     aYL, aT2L, aDmax, aT1max, aT2max, iC);
1093     }
1094     //
1095     if (iErr) {
1096       if ((iErr == 2) && !bMaxDist) {
1097         aXP = (aA + aB) * 0.5;
1098         DistPC(aXP, theC1, theCriteria, theProjPC, 
1099                aYP, aT2P, aDmax, aT1max, aT2max, iC);
1100       }
1101       return iErr;
1102     }
1103     //
1104     if ((aB - aA) < theEps) {
1105       break;
1106     }
1107   }// for (;;) {
1108   //
1109   return iErr;
1110 }
1111 //=======================================================================
1112 //function : DistPC
1113 //purpose  : 
1114 //=======================================================================
1115 Standard_Integer DistPC(const Standard_Real aT1, 
1116                         const Handle(Geom_Curve)& theC1,
1117                         const Standard_Real theCriteria,
1118                         GeomAPI_ProjectPointOnCurve& theProjPC, 
1119                         Standard_Real& aD, 
1120                         Standard_Real& aT2,
1121                         Standard_Real& aDmax,
1122                         Standard_Real& aT1max,
1123                         Standard_Real& aT2max,
1124                         const Standard_Integer iC)
1125 {
1126   Standard_Integer iErr;
1127   //
1128   iErr = DistPC(aT1, theC1, theCriteria, theProjPC, aD, aT2, iC);
1129   if (iErr == 1) {
1130     return iErr;
1131   }
1132   //
1133   if (iC*(aD - aDmax) > 0) {
1134     aDmax = aD;
1135     aT1max = aT1;
1136     aT2max = aT2;
1137   }
1138   //
1139   return iErr;
1140 }
1141 //=======================================================================
1142 //function : DistPC
1143 //purpose  : 
1144 //=======================================================================
1145 Standard_Integer DistPC(const Standard_Real aT1, 
1146                         const Handle(Geom_Curve)& theC1,
1147                         const Standard_Real theCriteria, 
1148                         GeomAPI_ProjectPointOnCurve& theProjPC,
1149                         Standard_Real& aD, 
1150                         Standard_Real& aT2,
1151                         const Standard_Integer iC) 
1152 {
1153   Standard_Integer iErr, aNbP2;
1154   gp_Pnt aP1;
1155   //
1156   iErr = 0;
1157   theC1->D0(aT1, aP1);
1158   //
1159   theProjPC.Perform(aP1);
1160   aNbP2 = theProjPC.NbPoints();
1161   if (!aNbP2) {
1162     iErr = 1;// the point from aC1 can not be projected on aC2
1163     return iErr;
1164   }
1165   //
1166   aD  = theProjPC.LowerDistance();
1167   aT2 = theProjPC.LowerDistanceParameter();
1168   if (iC*(aD - theCriteria) > 0) {
1169     iErr = 2;// the distance is too big or small
1170   }
1171   //
1172   return iErr;
1173 }
1174
1175 //=======================================================================
1176 //function : SplitRangeOnSegments
1177 //purpose  : 
1178 //=======================================================================
1179 Standard_Integer SplitRangeOnSegments(const Standard_Real aT1, 
1180                                       const Standard_Real aT2,
1181                                       const Standard_Real theResolution,
1182                                       const Standard_Integer theNbSeg,
1183                                       IntTools_SequenceOfRanges& theSegments)
1184 {
1185   Standard_Real aDiff = aT2 - aT1;
1186   if (aDiff < theResolution || theNbSeg == 1) {
1187     theSegments.Append(IntTools_Range(aT1, aT2));
1188     return 1;
1189   }
1190   //
1191   Standard_Real aDt, aT1x, aT2x, aSeg;
1192   Standard_Integer aNbSegments, i;
1193   //
1194   aNbSegments = theNbSeg;
1195   aDt = aDiff / aNbSegments;
1196   if (aDt < theResolution) {
1197     aSeg = aDiff / theResolution;
1198     aNbSegments = Standard_Integer(aSeg) + 1;
1199     aDt = aDiff / aNbSegments;
1200   }
1201   //
1202   aT1x = aT1;
1203   for (i = 1; i < aNbSegments; ++i) {
1204     aT2x = aT1x + aDt;
1205     //
1206     IntTools_Range aR(aT1x, aT2x);
1207     theSegments.Append(aR);
1208     //
1209     aT1x = aT2x;
1210   }
1211   //
1212   IntTools_Range aR(aT1x, aT2);
1213   theSegments.Append(aR);
1214   //
1215   return aNbSegments;
1216 }
1217
1218 //=======================================================================
1219 //function : BndBuildBox
1220 //purpose  : 
1221 //=======================================================================
1222 void BndBuildBox(const BRepAdaptor_Curve& theBAC,
1223                  const Standard_Real aT1,
1224                  const Standard_Real aT2,
1225                  const Standard_Real theTol,
1226                  Bnd_Box& theBox)
1227 {
1228   Bnd_Box aB;
1229   BndLib_Add3dCurve::Add(theBAC, aT1, aT2, theTol, aB);
1230   theBox = aB;
1231 }
1232
1233 //=======================================================================
1234 //function : PointBoxDistance
1235 //purpose  : 
1236 //=======================================================================
1237 Standard_Real PointBoxDistance(const Bnd_Box& aB,
1238                                const gp_Pnt& aP)
1239 {
1240   Standard_Real aPCoord[3];
1241   Standard_Real aBMinCoord[3], aBMaxCoord[3];
1242   Standard_Real aDist, aR1, aR2;
1243   Standard_Integer i;
1244   //
1245   aP.Coord(aPCoord[0], aPCoord[1], aPCoord[2]);
1246   aB.Get(aBMinCoord[0], aBMinCoord[1], aBMinCoord[2], 
1247          aBMaxCoord[0], aBMaxCoord[1], aBMaxCoord[2]);
1248   //
1249   aDist = 0.;
1250   for (i = 0; i < 3; ++i) {
1251     aR1 = aBMinCoord[i] - aPCoord[i];
1252     if (aR1 > 0.) {
1253       aDist += aR1*aR1;
1254       continue;
1255     }
1256     //
1257     aR2 = aPCoord[i] - aBMaxCoord[i];
1258     if (aR2 > 0.) {
1259       aDist += aR2*aR2;
1260     }
1261   }
1262   //
1263   aDist = Sqrt(aDist);
1264   return aDist;
1265 }
1266
1267 //=======================================================================
1268 //function : TypeToInteger
1269 //purpose  : 
1270 //=======================================================================
1271 Standard_Integer TypeToInteger(const GeomAbs_CurveType theCType)
1272 {
1273   Standard_Integer iRet;
1274   //
1275   switch(theCType) {
1276   case GeomAbs_Line:
1277     iRet=0;
1278     break;
1279   case GeomAbs_Hyperbola:
1280   case GeomAbs_Parabola:
1281     iRet=1;
1282     break;
1283   case GeomAbs_Circle:
1284   case GeomAbs_Ellipse:
1285     iRet=2;
1286     break;
1287   case GeomAbs_BezierCurve:
1288   case GeomAbs_BSplineCurve:
1289     iRet=3;
1290     break;
1291   default:
1292     iRet=4;
1293     break;
1294   }
1295   return iRet;
1296 }
1297
1298 //=======================================================================
1299 //function : ResolutionCoeff
1300 //purpose  : 
1301 //=======================================================================
1302 Standard_Real ResolutionCoeff(const BRepAdaptor_Curve& theBAC,
1303                               const IntTools_Range& theRange)
1304 {
1305   Standard_Real aResCoeff;
1306   //
1307   const Handle(Geom_Curve)& aCurve = theBAC.Curve().Curve();
1308   const GeomAbs_CurveType aCurveType = theBAC.GetType();
1309   //
1310   switch (aCurveType) {
1311   case GeomAbs_Circle :
1312     aResCoeff = 1. / (2 * Handle(Geom_Circle)::DownCast (aCurve)->Circ().Radius());
1313     break;
1314   case GeomAbs_Ellipse :
1315     aResCoeff =  1. / Handle(Geom_Ellipse)::DownCast (aCurve)->MajorRadius();
1316     break;
1317   case GeomAbs_Hyperbola :
1318   case GeomAbs_Parabola : 
1319   case GeomAbs_OffsetCurve :
1320   case GeomAbs_OtherCurve :{
1321     Standard_Real k, kMin, aDist, aDt, aT1, aT2, aT;
1322     Standard_Integer aNbP, i;
1323     gp_Pnt aP1, aP2;
1324     //
1325     aNbP = 30;
1326     theRange.Range(aT1, aT2);
1327     aDt = (aT2 - aT1) / aNbP;
1328     aT = aT1;
1329     kMin = 10.;
1330     //
1331     theBAC.D0(aT1, aP1);
1332     for (i = 1; i <= aNbP; ++i) {
1333       aT += aDt;
1334       theBAC.D0(aT, aP2);
1335       aDist = aP1.Distance(aP2);
1336       k = aDt / aDist;
1337       if (k < kMin) {
1338         kMin = k;
1339       }
1340       aP1 = aP2;
1341     }
1342     //
1343     aResCoeff = kMin;
1344     break;
1345   }
1346   default:
1347     aResCoeff = 0.;
1348     break;
1349   }
1350   //
1351   return aResCoeff;
1352 }
1353
1354 //=======================================================================
1355 //function : Resolution
1356 //purpose  : 
1357 //=======================================================================
1358 Standard_Real Resolution(const Handle(Geom_Curve)& theCurve,
1359                          const GeomAbs_CurveType theCurveType,
1360                          const Standard_Real theResCoeff,
1361                          const Standard_Real theR3D)
1362 {
1363   Standard_Real aRes;
1364   //
1365   switch (theCurveType) {
1366   case GeomAbs_Line :
1367     aRes = theR3D;
1368     break;
1369   case GeomAbs_Circle: {
1370     Standard_Real aDt = theResCoeff * theR3D;
1371     aRes = (aDt <= 1.) ? 2*ASin(aDt) : 2*M_PI;
1372     break;
1373   }
1374   case GeomAbs_BezierCurve:
1375     Handle(Geom_BezierCurve)::DownCast (theCurve)->Resolution(theR3D, aRes);
1376     break;
1377   case GeomAbs_BSplineCurve:
1378     Handle(Geom_BSplineCurve)::DownCast (theCurve)->Resolution(theR3D, aRes);
1379     break;
1380   default:
1381     aRes = theResCoeff * theR3D;
1382     break;
1383   }
1384   //
1385   return aRes;
1386 }
1387
1388 //=======================================================================
1389 //function : CurveDeflection
1390 //purpose  : 
1391 //=======================================================================
1392 Standard_Real CurveDeflection(const BRepAdaptor_Curve& theBAC,
1393                               const IntTools_Range& theRange)
1394 {
1395   Standard_Real aDt, aT, aT1, aT2, aDefl;
1396   Standard_Integer i, aNbP;
1397   gp_Vec aV1, aV2;
1398   gp_Pnt aP;
1399   //
1400   aDefl = 0;
1401   aNbP = 10;
1402   theRange.Range(aT1, aT2);
1403   aDt = (aT2 - aT1) / aNbP;
1404   aT = aT1;
1405   //
1406   theBAC.D1(aT1, aP, aV1);
1407   for (i = 1; i <= aNbP; ++i) {
1408     aT += aDt;
1409     theBAC.D1(aT, aP, aV2);
1410     if (aV1.Magnitude() > gp::Resolution() &&
1411         aV2.Magnitude() > gp::Resolution()) {
1412       gp_Dir aD1(aV1), aD2(aV2);
1413       aDefl += aD1.Angle(aD2);
1414     }
1415     aV1 = aV2;
1416   }
1417   //
1418   return aDefl;
1419 }
1420
1421 //=======================================================================
1422 //function : IsClosed
1423 //purpose  : 
1424 //=======================================================================
1425 Standard_Integer IsClosed(const Handle(Geom_Curve)& theCurve,
1426                           const Standard_Real aT1,
1427                           const Standard_Real aT2,
1428                           const Standard_Real theTol,
1429                           const Standard_Real theRes)
1430 {
1431   Standard_Boolean bClosed;
1432   Standard_Real aD;
1433   gp_Pnt aP1, aP2;
1434   //
1435   bClosed = Standard_False;
1436   if (Abs(aT1 - aT2) < theRes) {
1437     return bClosed;
1438   }
1439   //
1440   theCurve->D0(aT1, aP1);
1441   theCurve->D0(aT2, aP2);
1442   //
1443   aD = aP1.Distance(aP2);
1444   bClosed = aD < theTol;
1445   //
1446   return bClosed;
1447 }