1 // Created by: Eugeny MALTCHIKOV
2 // Copyright (c) 2013-2014 OPEN CASCADE SAS
4 // This file is part of Open CASCADE Technology software library.
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.
12 // Alternatively, this file may be used under the terms of Open CASCADE
13 // commercial license or contractual agreement.
15 #include <IntTools_EdgeEdge.ixx>
22 #include <TopoDS_Iterator.hxx>
24 #include <Bnd_Box.hxx>
25 #include <BndLib_Add3dCurve.hxx>
27 #include <Geom_Circle.hxx>
28 #include <Geom_Ellipse.hxx>
29 #include <Geom_BezierCurve.hxx>
30 #include <Geom_BSplineCurve.hxx>
32 #include <GeomAPI_ProjectPointOnCurve.hxx>
34 #include <BRep_Tool.hxx>
35 #include <BRepAdaptor_Curve.hxx>
37 #include <IntTools_CommonPrt.hxx>
39 #include <BOPCol_MapOfInteger.hxx>
42 void BndBuildBox(const BRepAdaptor_Curve& theBAC,
43 const Standard_Real aT1,
44 const Standard_Real aT2,
45 const Standard_Real theTol,
48 Standard_Real PointBoxDistance(const Bnd_Box& aB,
51 Standard_Integer SplitRangeOnSegments(const Standard_Real aT1,
52 const Standard_Real aT2,
53 const Standard_Real theResolution,
54 const Standard_Integer theNbSeg,
55 IntTools_SequenceOfRanges& theSegments);
57 Standard_Integer DistPC(const Standard_Real aT1,
58 const Handle(Geom_Curve)& theC1,
59 const Standard_Real theCriteria,
60 GeomAPI_ProjectPointOnCurve& theProjector,
63 const Standard_Integer iC = 1);
65 Standard_Integer DistPC(const Standard_Real aT1,
66 const Handle(Geom_Curve)& theC1,
67 const Standard_Real theCriteria,
68 GeomAPI_ProjectPointOnCurve& theProjector,
72 Standard_Real& aT1max,
73 Standard_Real& aT2max,
74 const Standard_Integer iC = 1);
76 Standard_Integer FindDistPC(const Standard_Real aT1A,
77 const Standard_Real aT1B,
78 const Handle(Geom_Curve)& theC1,
79 const Standard_Real theCriteria,
80 const Standard_Real theEps,
81 GeomAPI_ProjectPointOnCurve& theProjector,
83 Standard_Real& aT1max,
84 Standard_Real& aT2max,
85 const Standard_Boolean bMaxDist = Standard_True);
87 Standard_Real ResolutionCoeff(const BRepAdaptor_Curve& theBAC,
88 const IntTools_Range& theRange);
90 Standard_Real Resolution(const Handle(Geom_Curve)& theCurve,
91 const GeomAbs_CurveType theCurveType,
92 const Standard_Real theResCoeff,
93 const Standard_Real theR3D);
95 Standard_Real CurveDeflection(const BRepAdaptor_Curve& theBAC,
96 const IntTools_Range& theRange);
98 Standard_Integer IsClosed(const Handle(Geom_Curve)& theCurve,
99 const Standard_Real aT1,
100 const Standard_Real aT2,
101 const Standard_Real theTol,
102 const Standard_Real theRes);
104 Standard_Integer TypeToInteger(const GeomAbs_CurveType theCType);
106 //=======================================================================
109 //=======================================================================
110 void IntTools_EdgeEdge::Prepare()
112 GeomAbs_CurveType aCT1, aCT2;
113 Standard_Integer iCT1, iCT2;
115 myCurve1.Initialize(myEdge1);
116 myCurve2.Initialize(myEdge2);
118 if (myRange1.First() == 0. && myRange1.Last() == 0.) {
119 myRange1.SetFirst(myCurve1.FirstParameter());
120 myRange1.SetLast (myCurve1.LastParameter());
123 if (myRange2.First() == 0. && myRange2.Last() == 0.) {
124 myRange2.SetFirst(myCurve2.FirstParameter());
125 myRange2.SetLast (myCurve2.LastParameter());
128 aCT1 = myCurve1.GetType();
129 aCT2 = myCurve2.GetType();
131 iCT1 = TypeToInteger(aCT1);
132 iCT2 = TypeToInteger(aCT2);
137 Standard_Real aC1, aC2;
139 aC2 = CurveDeflection(myCurve2, myRange2);
140 aC1 = (aC2 > Precision::Confusion()) ?
141 CurveDeflection(myCurve1, myRange1) : 1.;
150 TopoDS_Edge tmpE = myEdge1;
154 BRepAdaptor_Curve tmpC = myCurve1;
158 IntTools_Range tmpR = myRange1;
162 mySwap = Standard_True;
165 myTol1 = myCurve1.Tolerance();
166 myTol2 = myCurve2.Tolerance();
167 myTol = myTol1 + myTol2;
169 if (iCT1 != 0 || iCT2 != 0) {
170 Standard_Real f, l, aTM;
172 myGeom1 = BRep_Tool::Curve(myEdge1, f, l);
173 myGeom2 = BRep_Tool::Curve(myEdge2, f, l);
175 myResCoeff1 = ResolutionCoeff(myCurve1, myRange1);
176 myResCoeff2 = ResolutionCoeff(myCurve2, myRange2);
178 myRes1 = Resolution(myCurve1.Curve().Curve(), myCurve1.GetType(), myResCoeff1, myTol1);
179 myRes2 = Resolution(myCurve2.Curve().Curve(), myCurve2.GetType(), myResCoeff2, myTol2);
182 aTM = Max(fabs(myRange1.First()), fabs(myRange1.Last()));
184 myPTol1 = 5.e-16 * aTM;
188 aTM = Max(fabs(myRange2.First()), fabs(myRange2.Last()));
190 myPTol2 = 5.e-16 * aTM;
195 //=======================================================================
198 //=======================================================================
199 void IntTools_EdgeEdge::Perform()
210 //3.1. Check Line/Line case
211 if (myCurve1.GetType() == GeomAbs_Line &&
212 myCurve2.GetType() == GeomAbs_Line) {
217 IntTools_SequenceOfRanges aRanges1, aRanges2;
219 //3.2. Find ranges containig solutions
220 Standard_Boolean bSplit2;
221 FindSolutions(aRanges1, aRanges2, bSplit2);
223 //4. Merge solutions and save common parts
224 MergeSolutions(aRanges1, aRanges2, bSplit2);
227 //=======================================================================
228 //function : FindSolutions
230 //=======================================================================
231 void IntTools_EdgeEdge::FindSolutions(IntTools_SequenceOfRanges& theRanges1,
232 IntTools_SequenceOfRanges& theRanges2,
233 Standard_Boolean& bSplit2)
235 Standard_Boolean bIsClosed2;
236 Standard_Real aT11, aT12, aT21, aT22;
239 bSplit2 = Standard_False;
240 myRange1.Range(aT11, aT12);
241 myRange2.Range(aT21, aT22);
243 bIsClosed2 = IsClosed(myGeom2, aT21, aT22, myTol2, myRes2);
247 BndBuildBox(myCurve1, aT11, aT12, myTol1, aB1);
249 gp_Pnt aP = myGeom2->Value(aT21);
250 bIsClosed2 = !aB1.IsOut(aP);
254 BndBuildBox(myCurve2, aT21, aT22, myTol2, aB2);
255 FindSolutions(myRange1, myRange2, aB2, theRanges1, theRanges2);
259 if (!CheckCoincidence(aT11, aT12, aT21, aT22, myTol, myRes1)) {
260 theRanges1.Append(myRange1);
261 theRanges2.Append(myRange2);
265 Standard_Integer i, j, aNb1, aNb2;
266 IntTools_SequenceOfRanges aSegments1, aSegments2;
268 aNb1 = IsClosed(myGeom1, aT11, aT12, myTol1, myRes1) ? 2 : 1;
271 aNb1 = SplitRangeOnSegments(aT11, aT12, myRes1, aNb1, aSegments1);
272 aNb2 = SplitRangeOnSegments(aT21, aT22, myRes2, aNb2, aSegments2);
274 for (i = 1; i <= aNb1; ++i) {
275 const IntTools_Range& aR1 = aSegments1(i);
276 for (j = 1; j <= aNb2; ++j) {
277 const IntTools_Range& aR2 = aSegments2(j);
278 BndBuildBox(myCurve2, aR2.First(), aR2.Last(), myTol2, aB2);
279 FindSolutions(aR1, aR2, aB2, theRanges1, theRanges2);
286 //=======================================================================
287 //function : FindSolutions
289 //=======================================================================
290 void IntTools_EdgeEdge::FindSolutions(const IntTools_Range& theR1,
291 const IntTools_Range& theR2,
292 const Bnd_Box& theBox2,
293 IntTools_SequenceOfRanges& theRanges1,
294 IntTools_SequenceOfRanges& theRanges2)
296 Standard_Boolean bOut, bStop, bThin;
297 Standard_Real aT11, aT12, aT21, aT22;
298 Standard_Real aTB11, aTB12, aTB21, aTB22;
299 Standard_Real aSmallStep1, aSmallStep2;
300 Standard_Integer iCom;
303 theR1.Range(aT11, aT12);
304 theR2.Range(aT21, aT22);
308 bThin = Standard_False;
309 bStop = Standard_False;
318 //1. Build box for first edge and find parameters
319 // of the second one in that box
320 BndBuildBox(myCurve1, aT11, aT12, myTol1, aB1);
321 bOut = aB1.IsOut(aB2);
326 bThin = ((aT12 - aT11) < myRes1) ||
327 (aB1.IsXThin(myTol) && aB1.IsYThin(myTol) && aB1.IsZThin(myTol));
329 bOut = !FindParameters(myCurve2, aTB21, aTB22, myRes2, myPTol2,
330 myResCoeff2, aB1, aT21, aT22);
335 //2. Build box for second edge and find parameters
336 // of the first one in that box
337 BndBuildBox(myCurve2, aT21, aT22, myTol2, aB2);
338 bOut = aB1.IsOut(aB2);
343 bThin = ((aT22 - aT21) < myRes2) ||
344 (aB2.IsXThin(myTol) && aB2.IsYThin(myTol) && aB2.IsZThin(myTol));
346 bOut = !FindParameters(myCurve1, aTB11, aTB12, myRes1, myPTol1,
347 myResCoeff1, aB2, aT11, aT12);
353 //3. Check if it makes sense to continue
354 aSmallStep1 = (aTB12 - aTB11) / 250.;
355 aSmallStep2 = (aTB22 - aTB21) / 250.;
357 if (aSmallStep1 < myRes1) {
358 aSmallStep1 = myRes1;
360 if (aSmallStep2 < myRes2) {
361 aSmallStep2 = myRes2;
364 if (((aT11 - aTB11) < aSmallStep1) && ((aTB12 - aT12) < aSmallStep1) &&
365 ((aT21 - aTB21) < aSmallStep2) && ((aTB22 - aT22) < aSmallStep2)) {
366 bStop = Standard_True;
377 //check curves for coincidence on the ranges
378 iCom = CheckCoincidence(aT11, aT12, aT21, aT22, myTol, myRes1);
380 bThin = Standard_True;
386 //check intermediate points
387 Standard_Boolean bSol;
390 GeomAPI_ProjectPointOnCurve aProjPC;
392 aT1 = (aT11 + aT12) * .5;
393 myGeom1->D0(aT1, aP1);
395 aProjPC.Init(myGeom2, aT21, aT22);
396 aProjPC.Perform(aP1);
398 if (aProjPC.NbPoints()) {
399 bSol = aProjPC.LowerDistance() <= myTol;
405 aT2 = (aT21 + aT22) * .5;
406 myGeom2->D0(aT2, aP2);
408 bSol = aP1.IsEqual(aP2, myTol);
416 IntTools_Range aR1(aT11, aT12), aR2(aT21, aT22);
418 theRanges1.Append(aR1);
419 theRanges2.Append(aR2);
423 if (!IsIntersection(aT11, aT12, aT21, aT22)) {
427 //split ranges on segments and repeat
428 Standard_Integer i, aNb1;
429 IntTools_SequenceOfRanges aSegments1;
431 IntTools_Range aR2(aT21, aT22);
432 BndBuildBox(myCurve2, aT21, aT22, myTol2, aB2);
434 aNb1 = SplitRangeOnSegments(aT11, aT12, myRes1, 3, aSegments1);
435 for (i = 1; i <= aNb1; ++i) {
436 const IntTools_Range& aR1 = aSegments1(i);
437 FindSolutions(aR1, aR2, aB2, theRanges1, theRanges2);
441 //=======================================================================
442 //function : FindParameters
444 //=======================================================================
445 Standard_Boolean IntTools_EdgeEdge::FindParameters(const BRepAdaptor_Curve& theBAC,
446 const Standard_Real aT1,
447 const Standard_Real aT2,
448 const Standard_Real theRes,
449 const Standard_Real thePTol,
450 const Standard_Real theResCoeff,
451 const Bnd_Box& theCBox,
455 Standard_Boolean bRet;
456 Standard_Integer aC, i, k;
457 Standard_Real aCf, aDiff, aDt, aT, aTB, aTOut, aTIn;
458 Standard_Real aDist, aDistP, aDistTol, aTol;
462 bRet = Standard_False;
463 aCf = 0.6180339887498948482045868343656;// =0.5*(1.+sqrt(5.))/2.;
465 aTol = theBAC.Tolerance();
467 aDistTol = Precision::PConfusion();
471 const Handle(Geom_Curve)& aCurve = theBAC.Curve().Curve();
472 const GeomAbs_CurveType aCurveType = theBAC.GetType();
474 for (i = 0; i < 2; ++i) {
475 aTB = !i ? aT1 : aT2;
476 aT = !i ? aT2 : aTB1;
478 bRet = Standard_False;
480 //looking for the point on the edge which is in the box;
481 while (aC*(aT-aTB) >= 0) {
483 aDist = PointBoxDistance(theCBox, aP);
485 if (fabs(aDist - aDistP) < aDistTol) {
486 aDt = Resolution(aCurve, aCurveType, theResCoeff, (++k)*aDist);
489 aDt = Resolution(aCurve, aCurveType, theResCoeff, aDist);
493 bRet = Standard_True;
501 //edge is out of the box;
512 //one point IN, one point OUT; looking for the bounding point;
514 aTOut = aTB - aC*aDt;
515 aDiff = aTIn - aTOut;
516 while (fabs(aDiff) > thePTol) {
517 aTB = aTOut + aDiff*aCf;
519 if (aCBx.IsOut(aP)) {
524 aDiff = aTIn - aTOut;
536 //=======================================================================
537 //function : MergeSolutions
539 //=======================================================================
540 void IntTools_EdgeEdge::MergeSolutions(const IntTools_SequenceOfRanges& theRanges1,
541 const IntTools_SequenceOfRanges& theRanges2,
542 const Standard_Boolean bSplit2)
544 Standard_Integer aNbCP = theRanges1.Length();
549 IntTools_Range aRi1, aRi2, aRj1, aRj2;
550 Standard_Boolean bCond;
551 Standard_Integer i, j;
552 TopAbs_ShapeEnum aType;
553 Standard_Real aT11, aT12, aT21, aT22;
554 Standard_Real aTi11, aTi12, aTi21, aTi22;
555 Standard_Real aTj11, aTj12, aTj21, aTj22;
556 Standard_Real aRes1, aRes2, dTR1, dTR2;
557 BOPCol_MapOfInteger aMI;
559 aRes1 = Resolution(myCurve1.Curve().Curve(),
560 myCurve1.GetType(), myResCoeff1, myTol);
561 aRes2 = Resolution(myCurve2.Curve().Curve(),
562 myCurve2.GetType(), myResCoeff2, myTol);
564 myRange1.Range(aT11, aT12);
565 myRange2.Range(aT21, aT22);
568 aType = TopAbs_VERTEX;
570 for (i = 1; i <= aNbCP;) {
571 if (aMI.Contains(i)) {
576 aRi1 = theRanges1(i);
577 aRi2 = theRanges2(i);
579 aRi1.Range(aTi11, aTi12);
580 aRi2.Range(aTi21, aTi22);
584 for (j = i+1; j <= aNbCP; ++j) {
585 if (aMI.Contains(j)) {
589 aRj1 = theRanges1(j);
590 aRj2 = theRanges2(j);
592 aRj1.Range(aTj11, aTj12);
593 aRj2.Range(aTj21, aTj22);
595 bCond = (fabs(aTi12 - aTj11) < dTR1) ||
596 (bSplit2 && (fabs(aTj12 - aTi11) < dTR1));
597 if (bCond && bSplit2) {
598 bCond = (fabs((Max(aTi22, aTj22) - Min(aTi21, aTj21)) -
599 ((aTi22 - aTi21) + (aTj22 - aTj21))) < dTR2);
603 aTi11 = Min(aTi11, aTj11);
604 aTi12 = Max(aTi12, aTj12);
605 aTi21 = Min(aTi21, aTj21);
606 aTi22 = Max(aTi22, aTj22);
615 if (((fabs(aT11 - aTi11) < myRes1) && (fabs(aT12 - aTi12) < myRes1)) ||
616 ((fabs(aT21 - aTi21) < myRes2) && (fabs(aT22 - aTi22) < myRes2))) {
618 myCommonParts.Clear();
621 AddSolution(aTi11, aTi12, aTi21, aTi22, aType);
622 if (aType == TopAbs_EDGE) {
632 //=======================================================================
633 //function : AddSolution
635 //=======================================================================
636 void IntTools_EdgeEdge::AddSolution(const Standard_Real aT11,
637 const Standard_Real aT12,
638 const Standard_Real aT21,
639 const Standard_Real aT22,
640 const TopAbs_ShapeEnum theType)
642 IntTools_CommonPrt aCPart;
644 aCPart.SetType(theType);
646 aCPart.SetEdge1(myEdge1);
647 aCPart.SetEdge2(myEdge2);
648 aCPart.SetRange1(aT11, aT12);
649 aCPart.AppendRange2(aT21, aT22);
651 aCPart.SetEdge1(myEdge2);
652 aCPart.SetEdge2(myEdge1);
653 aCPart.SetRange1(aT21, aT22);
654 aCPart.AppendRange2(aT11, aT12);
657 if (theType == TopAbs_VERTEX) {
658 Standard_Real aT1, aT2;
660 FindBestSolution(aT11, aT12, aT21, aT22, aT1, aT2);
663 aCPart.SetVertexParameter1(aT1);
664 aCPart.SetVertexParameter2(aT2);
666 aCPart.SetVertexParameter1(aT2);
667 aCPart.SetVertexParameter2(aT1);
670 myCommonParts.Append(aCPart);
673 //=======================================================================
674 //function : FindBestSolution
676 //=======================================================================
677 void IntTools_EdgeEdge::FindBestSolution(const Standard_Real aT11,
678 const Standard_Real aT12,
679 const Standard_Real aT21,
680 const Standard_Real aT22,
684 Standard_Integer i, aNbS, iErr;
685 Standard_Real aDMin, aD, aRes1, aSolCriteria, aTouchCriteria;
686 Standard_Real aT1A, aT1B, aT1Min, aT2Min;
687 Standard_Real aT1Im, aT2Im, aT1Touch;
688 GeomAPI_ProjectPointOnCurve aProjPC;
689 IntTools_SequenceOfRanges aRanges;
690 Standard_Boolean bTouch;
692 aDMin = Precision::Infinite();
693 aSolCriteria = 5.e-16;
694 aTouchCriteria = 5.e-13;
695 bTouch = Standard_False;
698 aRes1 = Resolution(myCurve1.Curve().Curve(),
699 myCurve1.GetType(), myResCoeff1, myTol);
701 aNbS = SplitRangeOnSegments(aT11, aT12, 3*aRes1, aNbS, aRanges);
703 aProjPC.Init(myGeom2, aT21, aT22);
705 aT1 = (aT11 + aT12) * 0.5;
706 iErr = DistPC(aT1, myGeom1, aSolCriteria, aProjPC, aD, aT2, -1);
708 aT2 = (aT21 + aT22) * 0.5;
714 for (i = 1; i <= aNbS; ++i) {
715 const IntTools_Range& aR1 = aRanges(i);
716 aR1.Range(aT1A, aT1B);
719 iErr = FindDistPC(aT1A, aT1B, myGeom1, aSolCriteria, myPTol1,
720 aProjPC, aD, aT1Min, aT2Min, Standard_False);
728 if (aD < aTouchCriteria) {
730 aT1A = (aT1Touch + aT1Min) * 0.5;
731 iErr = DistPC(aT1A, myGeom1, aTouchCriteria,
732 aProjPC, aD, aT2Min, -1);
733 if (aD > aTouchCriteria) {
741 bTouch = Standard_True;
748 //=======================================================================
749 //function : ComputeLineLine
751 //=======================================================================
752 void IntTools_EdgeEdge::ComputeLineLine()
754 Standard_Boolean IsParallel, IsCoincide;
755 Standard_Real aSin, aCos, aAng, aTol;
756 Standard_Real aT1, aT2, aT11, aT12, aT21, aT22;
760 IntTools_CommonPrt aCommonPrt;
762 IsParallel = Standard_False;
763 IsCoincide = Standard_False;
765 aL1 = myCurve1.Line();
766 aL2 = myCurve2.Line();
767 aD1 = aL1.Position().Direction();
768 aD2 = aL2.Position().Direction();
769 myRange1.Range(aT11, aT12);
770 myRange2.Range(aT21, aT22);
772 aCommonPrt.SetEdge1(myEdge1);
773 aCommonPrt.SetEdge2(myEdge2);
776 aAng = (aCos >= 0.) ? 2.*(1. - aCos) : 2.*(1. + aCos);
778 if(aAng <= Precision::Angular()) {
779 IsParallel = Standard_True;
780 if(aL1.SquareDistance(aL2.Location()) <= aTol) {
781 IsCoincide = Standard_True;
782 aP11 = ElCLib::Value(aT11, aL1);
783 aP12 = ElCLib::Value(aT12, aL1);
787 aP11 = ElCLib::Value(aT11, aL1);
788 aP12 = ElCLib::Value(aT12, aL1);
789 if(aL2.SquareDistance(aP11) <= aTol && aL2.SquareDistance(aP12) <= aTol) {
790 IsCoincide = Standard_True;
795 Standard_Real t21, t22;
797 t21 = ElCLib::Parameter(aL2, aP11);
798 t22 = ElCLib::Parameter(aL2, aP12);
799 if((t21 > aT22 && t22 > aT22) || (t21 < aT21 && t22 < aT21)) {
812 aCommonPrt.SetRange1(aT11, aT12);
813 aCommonPrt.SetAllNullFlag(Standard_True);
814 aCommonPrt.AppendRange2(t21, t22);
817 aCommonPrt.SetRange1(aT11, aT12 - (t22 - aT22));
818 aCommonPrt.AppendRange2(t21, aT22);
822 aCommonPrt.SetRange1(aT11 + (aT21 - t21), aT12);
823 aCommonPrt.AppendRange2(aT21, t22);
825 aCommonPrt.SetType(TopAbs_EDGE);
826 myCommonParts.Append(aCommonPrt);
835 TopoDS_Iterator aIt1, aIt2;
836 aIt1.Initialize(myEdge1);
837 for (; aIt1.More(); aIt1.Next()) {
838 const TopoDS_Shape& aV1 = aIt1.Value();
839 aIt2.Initialize(myEdge2);
840 for (; aIt2.More(); aIt2.Next()) {
841 const TopoDS_Shape& aV2 = aIt2.Value();
842 if (aV2.IsSame(aV1)) {
849 aSin = 1. - aCos*aCos;
850 gp_Pnt O1 = aL1.Location();
851 gp_Pnt O2 = aL2.Location();
852 gp_Vec O1O2 (O1, O2);
854 aT2 = (aD1.XYZ()*(O1O2.Dot(aD1))-(O1O2.XYZ())).Dot(aD2.XYZ());
857 if(aT2 < aT21 || aT2 > aT22) {
861 gp_Pnt aP2(ElCLib::Value(aT2, aL2));
862 aT1 = (gp_Vec(O1, aP2)).Dot(aD1);
864 if(aT1 < aT11 || aT1 > aT12) {
868 gp_Pnt aP1(ElCLib::Value(aT1, aL1));
869 Standard_Real aDist = aP1.SquareDistance(aP2);
875 aCommonPrt.SetRange1(aT1 - myTol, aT1 + myTol);
876 aCommonPrt.AppendRange2(aT2 - myTol, aT2 + myTol);
877 aCommonPrt.SetType(TopAbs_VERTEX);
878 aCommonPrt.SetVertexParameter1(aT1);
879 aCommonPrt.SetVertexParameter2(aT2);
880 myCommonParts.Append(aCommonPrt);
883 //=======================================================================
884 //function : IsIntersection
886 //=======================================================================
887 Standard_Boolean IntTools_EdgeEdge::IsIntersection(const Standard_Real aT11,
888 const Standard_Real aT12,
889 const Standard_Real aT21,
890 const Standard_Real aT22)
892 Standard_Boolean bRet;
893 gp_Pnt aP11, aP12, aP21, aP22;
894 gp_Vec aV11, aV12, aV21, aV22;
895 Standard_Real aD11_21, aD11_22, aD12_21, aD12_22, aCriteria, aCoef;
896 Standard_Boolean bSmall_11_21, bSmall_11_22, bSmall_12_21, bSmall_12_22;
898 bRet = Standard_True;
900 if (((aT12 - aT11) > aCoef*myRes1) && ((aT22 - aT21) > aCoef*myRes2)) {
903 Standard_Real aTRMin = Min((aT12 - aT11)/myRes1, (aT22 - aT21)/myRes2);
904 aCoef = aTRMin / 100.;
909 aCriteria = aCoef * myTol;
910 aCriteria *= aCriteria;
912 myGeom1->D1(aT11, aP11, aV11);
913 myGeom1->D1(aT12, aP12, aV12);
914 myGeom2->D1(aT21, aP21, aV21);
915 myGeom2->D1(aT22, aP22, aV22);
917 aD11_21 = aP11.SquareDistance(aP21);
918 aD11_22 = aP11.SquareDistance(aP22);
919 aD12_21 = aP12.SquareDistance(aP21);
920 aD12_22 = aP12.SquareDistance(aP22);
922 bSmall_11_21 = aD11_21 < aCriteria;
923 bSmall_11_22 = aD11_22 < aCriteria;
924 bSmall_12_21 = aD12_21 < aCriteria;
925 bSmall_12_22 = aD12_22 < aCriteria;
927 if ((bSmall_11_21 && bSmall_12_22) ||
928 (bSmall_11_22 && bSmall_12_21)) {
933 Standard_Real anAngleCriteria;
934 Standard_Real anAngle1, anAngle2;
936 anAngleCriteria = 5.e-3;
937 if (bSmall_11_21 && bSmall_12_22) {
938 anAngle1 = aV11.Angle(aV21);
939 anAngle2 = aV12.Angle(aV22);
941 anAngle1 = aV11.Angle(aV22);
942 anAngle2 = aV12.Angle(aV21);
945 if (((anAngle1 < anAngleCriteria) || ((M_PI - anAngle1) < anAngleCriteria)) ||
946 ((anAngle2 < anAngleCriteria) || ((M_PI - anAngle2) < anAngleCriteria))) {
947 GeomAPI_ProjectPointOnCurve aProjPC;
948 Standard_Integer iErr;
949 Standard_Real aD, aT1Min, aT2Min;
951 aD = Precision::Infinite();
952 aProjPC.Init(myGeom2, aT21, aT22);
953 iErr = FindDistPC(aT11, aT12, myGeom1, myTol, myRes1,
954 aProjPC, aD, aT1Min, aT2Min, Standard_False);
961 //=======================================================================
962 //function : CheckCoincidence
964 //=======================================================================
965 Standard_Integer IntTools_EdgeEdge::CheckCoincidence(const Standard_Real aT11,
966 const Standard_Real aT12,
967 const Standard_Real aT21,
968 const Standard_Real aT22,
969 const Standard_Real theCriteria,
970 const Standard_Real theCurveRes1)
972 Standard_Integer iErr, aNb, aNb1, i;
973 Standard_Real aT1A, aT1B, aT1max, aT2max, aDmax;
974 GeomAPI_ProjectPointOnCurve aProjPC;
975 IntTools_SequenceOfRanges aRanges;
979 aProjPC.Init(myGeom2, aT21, aT22);
981 // 1. Express evaluation
982 aNb = 10; // Number of intervals on the curve #1
983 aNb1 = SplitRangeOnSegments(aT11, aT12, theCurveRes1, aNb, aRanges);
984 for (i = 1; i < aNb1; ++i) {
985 const IntTools_Range& aR1 = aRanges(i);
986 aR1.Range(aT1A, aT1B);
988 iErr = DistPC(aT1B, myGeom1, theCriteria, aProjPC, aDmax, aT2max);
994 // if the ranges in aRanges are less than theCurveRes1,
995 // there is no need to do step 2 (deep evaluation)
1000 // 2. Deep evaluation
1001 for (i = 2; i < aNb1; ++i) {
1002 const IntTools_Range& aR1 = aRanges(i);
1003 aR1.Range(aT1A, aT1B);
1005 iErr = FindDistPC(aT1A, aT1B, myGeom1, theCriteria, theCurveRes1,
1006 aProjPC, aDmax, aT1max, aT2max);
1012 // iErr == 0 - the patches are coincided
1013 // iErr == 1 - a point from aC1 can not be projected on aC2
1014 // iErr == 2 - the distance is too big
1018 //=======================================================================
1019 //function : FindDistPC
1021 //=======================================================================
1022 Standard_Integer FindDistPC(const Standard_Real aT1A,
1023 const Standard_Real aT1B,
1024 const Handle(Geom_Curve)& theC1,
1025 const Standard_Real theCriteria,
1026 const Standard_Real theEps,
1027 GeomAPI_ProjectPointOnCurve& theProjPC,
1028 Standard_Real& aDmax,
1029 Standard_Real& aT1max,
1030 Standard_Real& aT2max,
1031 const Standard_Boolean bMaxDist)
1033 Standard_Integer iErr, iC;
1034 Standard_Real aGS, aXP, aA, aB, aXL, aYP, aYL, aT2P, aT2L;
1036 iC = bMaxDist ? 1 : -1;
1039 aGS = 0.6180339887498948482045868343656;// =0.5*(1.+sqrt(5.))-1.;
1044 iErr = DistPC(aA, theC1, theCriteria, theProjPC,
1045 aYP, aT2P, aDmax, aT1max, aT2max, iC);
1050 iErr = DistPC(aB, theC1, theCriteria, theProjPC,
1051 aYL, aT2L, aDmax, aT1max, aT2max, iC);
1056 aXP = aA + (aB - aA)*aGS;
1057 aXL = aB - (aB - aA)*aGS;
1059 iErr = DistPC(aXP, theC1, theCriteria, theProjPC,
1060 aYP, aT2P, aDmax, aT1max, aT2max, iC);
1065 iErr = DistPC(aXL, theC1, theCriteria, theProjPC,
1066 aYL, aT2L, aDmax, aT1max, aT2max, iC);
1072 if (iC*(aYP - aYL) > 0) {
1076 aXP = aA + (aB - aA)*aGS;
1077 iErr = DistPC(aXP, theC1, theCriteria, theProjPC,
1078 aYP, aT2P, aDmax, aT1max, aT2max, iC);
1084 aXL = aB - (aB - aA)*aGS;
1085 iErr = DistPC(aXL, theC1, theCriteria, theProjPC,
1086 aYL, aT2L, aDmax, aT1max, aT2max, iC);
1090 if ((iErr == 2) && !bMaxDist) {
1091 aXP = (aA + aB) * 0.5;
1092 DistPC(aXP, theC1, theCriteria, theProjPC,
1093 aYP, aT2P, aDmax, aT1max, aT2max, iC);
1098 if ((aB - aA) < theEps) {
1105 //=======================================================================
1108 //=======================================================================
1109 Standard_Integer DistPC(const Standard_Real aT1,
1110 const Handle(Geom_Curve)& theC1,
1111 const Standard_Real theCriteria,
1112 GeomAPI_ProjectPointOnCurve& theProjPC,
1115 Standard_Real& aDmax,
1116 Standard_Real& aT1max,
1117 Standard_Real& aT2max,
1118 const Standard_Integer iC)
1120 Standard_Integer iErr;
1122 iErr = DistPC(aT1, theC1, theCriteria, theProjPC, aD, aT2, iC);
1127 if (iC*(aD - aDmax) > 0) {
1135 //=======================================================================
1138 //=======================================================================
1139 Standard_Integer DistPC(const Standard_Real aT1,
1140 const Handle(Geom_Curve)& theC1,
1141 const Standard_Real theCriteria,
1142 GeomAPI_ProjectPointOnCurve& theProjPC,
1145 const Standard_Integer iC)
1147 Standard_Integer iErr, aNbP2;
1151 theC1->D0(aT1, aP1);
1153 theProjPC.Perform(aP1);
1154 aNbP2 = theProjPC.NbPoints();
1156 iErr = 1;// the point from aC1 can not be projected on aC2
1160 aD = theProjPC.LowerDistance();
1161 aT2 = theProjPC.LowerDistanceParameter();
1162 if (iC*(aD - theCriteria) > 0) {
1163 iErr = 2;// the distance is too big or small
1169 //=======================================================================
1170 //function : SplitRangeOnSegments
1172 //=======================================================================
1173 Standard_Integer SplitRangeOnSegments(const Standard_Real aT1,
1174 const Standard_Real aT2,
1175 const Standard_Real theResolution,
1176 const Standard_Integer theNbSeg,
1177 IntTools_SequenceOfRanges& theSegments)
1179 Standard_Real aDiff = aT2 - aT1;
1180 if (aDiff < theResolution || theNbSeg == 1) {
1181 theSegments.Append(IntTools_Range(aT1, aT2));
1185 Standard_Real aDt, aT1x, aT2x, aSeg;
1186 Standard_Integer aNbSegments, i;
1188 aNbSegments = theNbSeg;
1189 aDt = aDiff / aNbSegments;
1190 if (aDt < theResolution) {
1191 aSeg = aDiff / theResolution;
1192 aNbSegments = Standard_Integer(aSeg) + 1;
1193 aDt = aDiff / aNbSegments;
1197 for (i = 1; i < aNbSegments; ++i) {
1200 IntTools_Range aR(aT1x, aT2x);
1201 theSegments.Append(aR);
1206 IntTools_Range aR(aT1x, aT2);
1207 theSegments.Append(aR);
1212 //=======================================================================
1213 //function : BndBuildBox
1215 //=======================================================================
1216 void BndBuildBox(const BRepAdaptor_Curve& theBAC,
1217 const Standard_Real aT1,
1218 const Standard_Real aT2,
1219 const Standard_Real theTol,
1223 BndLib_Add3dCurve::Add(theBAC, aT1, aT2, theTol, aB);
1227 //=======================================================================
1228 //function : PointBoxDistance
1230 //=======================================================================
1231 Standard_Real PointBoxDistance(const Bnd_Box& aB,
1234 Standard_Real aPCoord[3];
1235 Standard_Real aBMinCoord[3], aBMaxCoord[3];
1236 Standard_Real aDist, aR1, aR2;
1239 aP.Coord(aPCoord[0], aPCoord[1], aPCoord[2]);
1240 aB.Get(aBMinCoord[0], aBMinCoord[1], aBMinCoord[2],
1241 aBMaxCoord[0], aBMaxCoord[1], aBMaxCoord[2]);
1244 for (i = 0; i < 3; ++i) {
1245 aR1 = aBMinCoord[i] - aPCoord[i];
1251 aR2 = aPCoord[i] - aBMaxCoord[i];
1257 aDist = Sqrt(aDist);
1261 //=======================================================================
1262 //function : TypeToInteger
1264 //=======================================================================
1265 Standard_Integer TypeToInteger(const GeomAbs_CurveType theCType)
1267 Standard_Integer iRet;
1273 case GeomAbs_Hyperbola:
1274 case GeomAbs_Parabola:
1277 case GeomAbs_Circle:
1278 case GeomAbs_Ellipse:
1281 case GeomAbs_BezierCurve:
1282 case GeomAbs_BSplineCurve:
1292 //=======================================================================
1293 //function : ResolutionCoeff
1295 //=======================================================================
1296 Standard_Real ResolutionCoeff(const BRepAdaptor_Curve& theBAC,
1297 const IntTools_Range& theRange)
1299 Standard_Real aResCoeff;
1301 const Handle(Geom_Curve)& aCurve = theBAC.Curve().Curve();
1302 const GeomAbs_CurveType aCurveType = theBAC.GetType();
1304 switch (aCurveType) {
1305 case GeomAbs_Circle :
1306 aResCoeff = 1. / (2 * Handle(Geom_Circle)::DownCast (aCurve)->Circ().Radius());
1308 case GeomAbs_Ellipse :
1309 aResCoeff = 1. / Handle(Geom_Ellipse)::DownCast (aCurve)->MajorRadius();
1311 case GeomAbs_Hyperbola :
1312 case GeomAbs_Parabola :
1313 case GeomAbs_OtherCurve :{
1314 Standard_Real k, kMin, aDist, aDt, aT1, aT2, aT;
1315 Standard_Integer aNbP, i;
1319 theRange.Range(aT1, aT2);
1320 aDt = (aT2 - aT1) / aNbP;
1324 theBAC.D0(aT1, aP1);
1325 for (i = 1; i <= aNbP; ++i) {
1328 aDist = aP1.Distance(aP2);
1347 //=======================================================================
1348 //function : Resolution
1350 //=======================================================================
1351 Standard_Real Resolution(const Handle(Geom_Curve)& theCurve,
1352 const GeomAbs_CurveType theCurveType,
1353 const Standard_Real theResCoeff,
1354 const Standard_Real theR3D)
1358 switch (theCurveType) {
1362 case GeomAbs_Circle: {
1363 Standard_Real aDt = theResCoeff * theR3D;
1364 aRes = (aDt <= 1.) ? 2*ASin(aDt) : 2*M_PI;
1367 case GeomAbs_BezierCurve:
1368 Handle(Geom_BezierCurve)::DownCast (theCurve)->Resolution(theR3D, aRes);
1370 case GeomAbs_BSplineCurve:
1371 Handle(Geom_BSplineCurve)::DownCast (theCurve)->Resolution(theR3D, aRes);
1374 aRes = theResCoeff * theR3D;
1381 //=======================================================================
1382 //function : CurveDeflection
1384 //=======================================================================
1385 Standard_Real CurveDeflection(const BRepAdaptor_Curve& theBAC,
1386 const IntTools_Range& theRange)
1388 Standard_Real aDt, aT, aT1, aT2, aDefl;
1389 Standard_Integer i, aNbP;
1395 theRange.Range(aT1, aT2);
1396 aDt = (aT2 - aT1) / aNbP;
1399 theBAC.D1(aT1, aP, aV1);
1400 for (i = 1; i <= aNbP; ++i) {
1402 theBAC.D1(aT, aP, aV2);
1403 if (aV1.Magnitude() > gp::Resolution() &&
1404 aV2.Magnitude() > gp::Resolution()) {
1405 gp_Dir aD1(aV1), aD2(aV2);
1406 aDefl += aD1.Angle(aD2);
1414 //=======================================================================
1415 //function : IsClosed
1417 //=======================================================================
1418 Standard_Integer IsClosed(const Handle(Geom_Curve)& theCurve,
1419 const Standard_Real aT1,
1420 const Standard_Real aT2,
1421 const Standard_Real theTol,
1422 const Standard_Real theRes)
1424 Standard_Boolean bClosed;
1428 bClosed = Standard_False;
1429 if (Abs(aT1 - aT2) < theRes) {
1433 theCurve->D0(aT1, aP1);
1434 theCurve->D0(aT2, aP2);
1436 aD = aP1.Distance(aP2);
1437 bClosed = aD < theTol;