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.
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>
22 #include <Geom_BezierCurve.hxx>
23 #include <Geom_BSplineCurve.hxx>
24 #include <Geom_Line.hxx>
25 #include <Geom_Circle.hxx>
26 #include <Geom_Curve.hxx>
27 #include <Geom_Ellipse.hxx>
28 #include <Geom_OffsetCurve.hxx>
29 #include <GeomAPI_ProjectPointOnCurve.hxx>
32 #include <IntTools_CommonPrt.hxx>
33 #include <IntTools_EdgeEdge.hxx>
34 #include <IntTools_Range.hxx>
35 #include <IntTools_Tools.hxx>
36 #include <TopoDS_Edge.hxx>
37 #include <TopoDS_Iterator.hxx>
40 void BndBuildBox(const BRepAdaptor_Curve& theBAC,
41 const Standard_Real aT1,
42 const Standard_Real aT2,
43 const Standard_Real theTol,
46 Standard_Real PointBoxDistance(const Bnd_Box& aB,
49 Standard_Integer SplitRangeOnSegments(const Standard_Real aT1,
50 const Standard_Real aT2,
51 const Standard_Real theResolution,
52 const Standard_Integer theNbSeg,
53 IntTools_SequenceOfRanges& theSegments);
55 Standard_Integer DistPC(const Standard_Real aT1,
56 const Handle(Geom_Curve)& theC1,
57 const Standard_Real theCriteria,
58 GeomAPI_ProjectPointOnCurve& theProjector,
61 const Standard_Integer iC = 1);
63 Standard_Integer DistPC(const Standard_Real aT1,
64 const Handle(Geom_Curve)& theC1,
65 const Standard_Real theCriteria,
66 GeomAPI_ProjectPointOnCurve& theProjector,
70 Standard_Real& aT1max,
71 Standard_Real& aT2max,
72 const Standard_Integer iC = 1);
74 Standard_Integer FindDistPC(const Standard_Real aT1A,
75 const Standard_Real aT1B,
76 const Handle(Geom_Curve)& theC1,
77 const Standard_Real theCriteria,
78 const Standard_Real theEps,
79 GeomAPI_ProjectPointOnCurve& theProjector,
81 Standard_Real& aT1max,
82 Standard_Real& aT2max,
83 const Standard_Boolean bMaxDist = Standard_True);
85 Standard_Real ResolutionCoeff(const BRepAdaptor_Curve& theBAC,
86 const IntTools_Range& theRange);
88 Standard_Real Resolution(const Handle(Geom_Curve)& theCurve,
89 const GeomAbs_CurveType theCurveType,
90 const Standard_Real theResCoeff,
91 const Standard_Real theR3D);
93 Standard_Real CurveDeflection(const BRepAdaptor_Curve& theBAC,
94 const IntTools_Range& theRange);
96 Standard_Boolean IsClosed(const Handle(Geom_Curve)& theCurve,
97 const Standard_Real aT1,
98 const Standard_Real aT2,
99 const Standard_Real theTol,
100 const Standard_Real theRes);
102 Standard_Integer TypeToInteger(const GeomAbs_CurveType theCType);
104 //=======================================================================
107 //=======================================================================
108 void IntTools_EdgeEdge::Prepare()
110 GeomAbs_CurveType aCT1, aCT2;
111 Standard_Integer iCT1, iCT2;
113 myCurve1.Initialize(myEdge1);
114 myCurve2.Initialize(myEdge2);
116 if (myRange1.First() == 0. && myRange1.Last() == 0.) {
117 myRange1.SetFirst(myCurve1.FirstParameter());
118 myRange1.SetLast (myCurve1.LastParameter());
121 if (myRange2.First() == 0. && myRange2.Last() == 0.) {
122 myRange2.SetFirst(myCurve2.FirstParameter());
123 myRange2.SetLast (myCurve2.LastParameter());
126 aCT1 = myCurve1.GetType();
127 aCT2 = myCurve2.GetType();
129 iCT1 = TypeToInteger(aCT1);
130 iCT2 = TypeToInteger(aCT2);
135 Standard_Real aC1, aC2;
137 aC2 = CurveDeflection(myCurve2, myRange2);
138 aC1 = (aC2 > Precision::Confusion()) ?
139 CurveDeflection(myCurve1, myRange1) : 1.;
148 TopoDS_Edge tmpE = myEdge1;
152 BRepAdaptor_Curve tmpC = myCurve1;
156 IntTools_Range tmpR = myRange1;
160 mySwap = Standard_True;
163 Standard_Real aTolAdd = myFuzzyValue / 2.;
164 myTol1 = myCurve1.Tolerance() + aTolAdd;
165 myTol2 = myCurve2.Tolerance() + aTolAdd;
166 myTol = myTol1 + myTol2;
168 if (iCT1 != 0 || iCT2 != 0) {
169 Standard_Real f, l, aTM;
171 myGeom1 = BRep_Tool::Curve(myEdge1, f, l);
172 myGeom2 = BRep_Tool::Curve(myEdge2, f, l);
174 myResCoeff1 = ResolutionCoeff(myCurve1, myRange1);
175 myResCoeff2 = ResolutionCoeff(myCurve2, myRange2);
177 myRes1 = Resolution(myCurve1.Curve().Curve(), myCurve1.GetType(), myResCoeff1, myTol1);
178 myRes2 = Resolution(myCurve2.Curve().Curve(), myCurve2.GetType(), myResCoeff2, myTol2);
181 aTM = Max(fabs(myRange1.First()), fabs(myRange1.Last()));
183 myPTol1 = 5.e-16 * aTM;
187 aTM = Max(fabs(myRange2.First()), fabs(myRange2.Last()));
189 myPTol2 = 5.e-16 * aTM;
194 //=======================================================================
197 //=======================================================================
198 void IntTools_EdgeEdge::Perform()
209 //3.1. Check Line/Line case
210 if (myCurve1.GetType() == GeomAbs_Line &&
211 myCurve2.GetType() == GeomAbs_Line) {
216 if (myQuickCoincidenceCheck) {
217 if (IsCoincident()) {
218 Standard_Real aT11, aT12, aT21, aT22;
220 myRange1.Range(aT11, aT12);
221 myRange2.Range(aT21, aT22);
222 AddSolution(aT11, aT12, aT21, aT22, TopAbs_EDGE);
227 IntTools_SequenceOfRanges aRanges1, aRanges2;
229 //3.2. Find ranges containig solutions
230 Standard_Boolean bSplit2;
231 FindSolutions(aRanges1, aRanges2, bSplit2);
233 //4. Merge solutions and save common parts
234 MergeSolutions(aRanges1, aRanges2, bSplit2);
237 //=======================================================================
238 //function : IsCoincident
240 //=======================================================================
241 Standard_Boolean IntTools_EdgeEdge::IsCoincident()
243 Standard_Integer i, iCnt, aNbSeg, aNbP2;
244 Standard_Real dT, aT1, aCoeff, aTresh, aD;
245 Standard_Real aT11, aT12, aT21, aT22;
246 GeomAPI_ProjectPointOnCurve aProjPC;
251 myRange1.Range(aT11, aT12);
252 myRange2.Range(aT21, aT22);
254 aProjPC.Init(myGeom2, aT21, aT22);
256 dT=(aT12-aT11)/aNbSeg;
259 for(i=0; i <= aNbSeg; ++i) {
261 myGeom1->D0(aT1, aP1);
263 aProjPC.Perform(aP1);
264 aNbP2=aProjPC.NbPoints();
269 aD=aProjPC.LowerDistance();
275 aCoeff=(Standard_Real)iCnt/((Standard_Real)aNbSeg+1);
276 return aCoeff > aTresh;
278 //=======================================================================
279 //function : FindSolutions
281 //=======================================================================
282 void IntTools_EdgeEdge::FindSolutions(IntTools_SequenceOfRanges& theRanges1,
283 IntTools_SequenceOfRanges& theRanges2,
284 Standard_Boolean& bSplit2)
286 Standard_Boolean bIsClosed2;
287 Standard_Real aT11, aT12, aT21, aT22;
290 bSplit2 = Standard_False;
291 myRange1.Range(aT11, aT12);
292 myRange2.Range(aT21, aT22);
294 bIsClosed2 = IsClosed(myGeom2, aT21, aT22, myTol2, myRes2);
298 BndBuildBox(myCurve1, aT11, aT12, myTol1, aB1);
300 gp_Pnt aP = myGeom2->Value(aT21);
301 bIsClosed2 = !aB1.IsOut(aP);
305 BndBuildBox(myCurve2, aT21, aT22, myTol2, aB2);
306 FindSolutions(myRange1, myRange2, aB2, theRanges1, theRanges2);
310 if (!CheckCoincidence(aT11, aT12, aT21, aT22, myTol, myRes1)) {
311 theRanges1.Append(myRange1);
312 theRanges2.Append(myRange2);
316 Standard_Integer i, j, aNb1, aNb2;
317 IntTools_SequenceOfRanges aSegments1, aSegments2;
319 aNb1 = IsClosed(myGeom1, aT11, aT12, myTol1, myRes1) ? 2 : 1;
322 aNb1 = SplitRangeOnSegments(aT11, aT12, myRes1, aNb1, aSegments1);
323 aNb2 = SplitRangeOnSegments(aT21, aT22, myRes2, aNb2, aSegments2);
325 for (i = 1; i <= aNb1; ++i) {
326 const IntTools_Range& aR1 = aSegments1(i);
327 for (j = 1; j <= aNb2; ++j) {
328 const IntTools_Range& aR2 = aSegments2(j);
329 BndBuildBox(myCurve2, aR2.First(), aR2.Last(), myTol2, aB2);
330 FindSolutions(aR1, aR2, aB2, theRanges1, theRanges2);
337 //=======================================================================
338 //function : FindSolutions
340 //=======================================================================
341 void IntTools_EdgeEdge::FindSolutions(const IntTools_Range& theR1,
342 const IntTools_Range& theR2,
343 const Bnd_Box& theBox2,
344 IntTools_SequenceOfRanges& theRanges1,
345 IntTools_SequenceOfRanges& theRanges2)
347 Standard_Boolean bOut, bStop, bThin;
348 Standard_Real aT11, aT12, aT21, aT22;
349 Standard_Real aTB11, aTB12, aTB21, aTB22;
350 Standard_Real aSmallStep1, aSmallStep2;
351 Standard_Integer iCom;
354 theR1.Range(aT11, aT12);
355 theR2.Range(aT21, aT22);
359 bThin = Standard_False;
360 bStop = Standard_False;
369 //1. Build box for first edge and find parameters
370 // of the second one in that box
371 BndBuildBox(myCurve1, aT11, aT12, myTol1, aB1);
372 bOut = aB1.IsOut(aB2);
377 bThin = ((aT12 - aT11) < myRes1) ||
378 (aB1.IsXThin(myTol) && aB1.IsYThin(myTol) && aB1.IsZThin(myTol));
380 bOut = !FindParameters(myCurve2, aTB21, aTB22, myTol2, myRes2, myPTol2,
381 myResCoeff2, aB1, aT21, aT22);
386 //2. Build box for second edge and find parameters
387 // of the first one in that box
388 BndBuildBox(myCurve2, aT21, aT22, myTol2, aB2);
389 bOut = aB1.IsOut(aB2);
394 bThin = ((aT22 - aT21) < myRes2) ||
395 (aB2.IsXThin(myTol) && aB2.IsYThin(myTol) && aB2.IsZThin(myTol));
397 bOut = !FindParameters(myCurve1, aTB11, aTB12, myTol1, myRes1, myPTol1,
398 myResCoeff1, aB2, aT11, aT12);
404 //3. Check if it makes sense to continue
405 aSmallStep1 = (aTB12 - aTB11) / 250.;
406 aSmallStep2 = (aTB22 - aTB21) / 250.;
408 if (aSmallStep1 < myRes1) {
409 aSmallStep1 = myRes1;
411 if (aSmallStep2 < myRes2) {
412 aSmallStep2 = myRes2;
415 if (((aT11 - aTB11) < aSmallStep1) && ((aTB12 - aT12) < aSmallStep1) &&
416 ((aT21 - aTB21) < aSmallStep2) && ((aTB22 - aT22) < aSmallStep2)) {
417 bStop = Standard_True;
428 //check curves for coincidence on the ranges
429 iCom = CheckCoincidence(aT11, aT12, aT21, aT22, myTol, myRes1);
431 bThin = Standard_True;
437 //check intermediate points
438 Standard_Boolean bSol;
441 GeomAPI_ProjectPointOnCurve aProjPC;
443 aT1 = (aT11 + aT12) * .5;
444 myGeom1->D0(aT1, aP1);
446 aProjPC.Init(myGeom2, aT21, aT22);
447 aProjPC.Perform(aP1);
449 if (aProjPC.NbPoints()) {
450 bSol = aProjPC.LowerDistance() <= myTol;
456 aT2 = (aT21 + aT22) * .5;
457 myGeom2->D0(aT2, aP2);
459 bSol = aP1.IsEqual(aP2, myTol);
467 IntTools_Range aR1(aT11, aT12), aR2(aT21, aT22);
469 theRanges1.Append(aR1);
470 theRanges2.Append(aR2);
474 if (!IsIntersection(aT11, aT12, aT21, aT22)) {
478 //split ranges on segments and repeat
479 Standard_Integer i, aNb1;
480 IntTools_SequenceOfRanges aSegments1;
482 IntTools_Range aR2(aT21, aT22);
483 BndBuildBox(myCurve2, aT21, aT22, myTol2, aB2);
485 aNb1 = SplitRangeOnSegments(aT11, aT12, myRes1, 3, aSegments1);
486 for (i = 1; i <= aNb1; ++i) {
487 const IntTools_Range& aR1 = aSegments1(i);
488 FindSolutions(aR1, aR2, aB2, theRanges1, theRanges2);
492 //=======================================================================
493 //function : FindParameters
495 //=======================================================================
496 Standard_Boolean IntTools_EdgeEdge::FindParameters(const BRepAdaptor_Curve& theBAC,
497 const Standard_Real aT1,
498 const Standard_Real aT2,
499 const Standard_Real theTol,
500 const Standard_Real theRes,
501 const Standard_Real thePTol,
502 const Standard_Real theResCoeff,
503 const Bnd_Box& theCBox,
507 Standard_Boolean bRet;
508 Standard_Integer aC, i;
509 Standard_Real aCf, aDiff, aDt, aT, aTB, aTOut, aTIn;
510 Standard_Real aDist, aDistP;
514 bRet = Standard_False;
515 aCf = 0.6180339887498948482045868343656;// =0.5*(1.+sqrt(5.))/2.;
517 aCBx.SetGap(aCBx.GetGap() + theTol);
519 const Handle(Geom_Curve)& aCurve = theBAC.Curve().Curve();
520 const GeomAbs_CurveType aCurveType = theBAC.GetType();
521 Standard_Real aMaxDt = (aT2 - aT1) * 0.01;
523 for (i = 0; i < 2; ++i) {
524 aTB = !i ? aT1 : aT2;
525 aT = !i ? aT2 : aTB1;
529 bRet = Standard_False;
531 //looking for the point on the edge which is in the box;
532 while (aC*(aT-aTB) >= 0) {
534 aDist = PointBoxDistance(theCBox, aP);
535 if (aDist > theTol) {
537 Standard_Boolean toGrow = Standard_False;
538 if (Abs(aDistP - aDist) / aDistP < 0.1) {
539 aDt = Resolution(aCurve, aCurveType, theResCoeff, k*aDist);
542 toGrow = Standard_True;
548 aDt = Resolution(aCurve, aCurveType, theResCoeff, aDist);
553 bRet = Standard_True;
561 //edge is out of the box;
572 //one point IN, one point OUT; looking for the bounding point;
574 aTOut = aTB - aC*aDt;
575 aDiff = aTIn - aTOut;
576 while (fabs(aDiff) > thePTol) {
577 aTB = aTOut + aDiff*aCf;
579 if (aCBx.IsOut(aP)) {
584 aDiff = aTIn - aTOut;
596 //=======================================================================
597 //function : MergeSolutions
599 //=======================================================================
600 void IntTools_EdgeEdge::MergeSolutions(const IntTools_SequenceOfRanges& theRanges1,
601 const IntTools_SequenceOfRanges& theRanges2,
602 const Standard_Boolean bSplit2)
604 Standard_Integer aNbCP = theRanges1.Length();
609 IntTools_Range aRi1, aRi2, aRj1, aRj2;
610 Standard_Boolean bCond;
611 Standard_Integer i, j;
612 TopAbs_ShapeEnum aType;
613 Standard_Real aT11, aT12, aT21, aT22;
614 Standard_Real aTi11, aTi12, aTi21, aTi22;
615 Standard_Real aTj11, aTj12, aTj21, aTj22;
616 Standard_Real aRes1, aRes2, dTR1, dTR2;
617 BOPCol_MapOfInteger aMI;
619 aRes1 = Resolution(myCurve1.Curve().Curve(),
620 myCurve1.GetType(), myResCoeff1, myTol);
621 aRes2 = Resolution(myCurve2.Curve().Curve(),
622 myCurve2.GetType(), myResCoeff2, myTol);
624 myRange1.Range(aT11, aT12);
625 myRange2.Range(aT21, aT22);
628 aType = TopAbs_VERTEX;
630 for (i = 1; i <= aNbCP;) {
631 if (aMI.Contains(i)) {
636 aRi1 = theRanges1(i);
637 aRi2 = theRanges2(i);
639 aRi1.Range(aTi11, aTi12);
640 aRi2.Range(aTi21, aTi22);
644 for (j = i+1; j <= aNbCP; ++j) {
645 if (aMI.Contains(j)) {
649 aRj1 = theRanges1(j);
650 aRj2 = theRanges2(j);
652 aRj1.Range(aTj11, aTj12);
653 aRj2.Range(aTj21, aTj22);
655 bCond = (fabs(aTi12 - aTj11) < dTR1) ||
656 (bSplit2 && (fabs(aTj12 - aTi11) < dTR1));
657 if (bCond && bSplit2) {
658 bCond = (fabs((Max(aTi22, aTj22) - Min(aTi21, aTj21)) -
659 ((aTi22 - aTi21) + (aTj22 - aTj21))) < dTR2);
663 aTi11 = Min(aTi11, aTj11);
664 aTi12 = Max(aTi12, aTj12);
665 aTi21 = Min(aTi21, aTj21);
666 aTi22 = Max(aTi22, aTj22);
675 if (((fabs(aT11 - aTi11) < myRes1) && (fabs(aT12 - aTi12) < myRes1)) ||
676 ((fabs(aT21 - aTi21) < myRes2) && (fabs(aT22 - aTi22) < myRes2))) {
678 myCommonParts.Clear();
681 AddSolution(aTi11, aTi12, aTi21, aTi22, aType);
682 if (aType == TopAbs_EDGE) {
692 //=======================================================================
693 //function : AddSolution
695 //=======================================================================
696 void IntTools_EdgeEdge::AddSolution(const Standard_Real aT11,
697 const Standard_Real aT12,
698 const Standard_Real aT21,
699 const Standard_Real aT22,
700 const TopAbs_ShapeEnum theType)
702 IntTools_CommonPrt aCPart;
704 aCPart.SetType(theType);
706 aCPart.SetEdge1(myEdge1);
707 aCPart.SetEdge2(myEdge2);
708 aCPart.SetRange1(aT11, aT12);
709 aCPart.AppendRange2(aT21, aT22);
711 aCPart.SetEdge1(myEdge2);
712 aCPart.SetEdge2(myEdge1);
713 aCPart.SetRange1(aT21, aT22);
714 aCPart.AppendRange2(aT11, aT12);
717 if (theType == TopAbs_VERTEX) {
718 Standard_Real aT1, aT2;
720 FindBestSolution(aT11, aT12, aT21, aT22, aT1, aT2);
723 aCPart.SetVertexParameter1(aT1);
724 aCPart.SetVertexParameter2(aT2);
726 aCPart.SetVertexParameter1(aT2);
727 aCPart.SetVertexParameter2(aT1);
730 myCommonParts.Append(aCPart);
733 //=======================================================================
734 //function : FindBestSolution
736 //=======================================================================
737 void IntTools_EdgeEdge::FindBestSolution(const Standard_Real aT11,
738 const Standard_Real aT12,
739 const Standard_Real aT21,
740 const Standard_Real aT22,
744 Standard_Integer i, aNbS, iErr;
745 Standard_Real aDMin, aD, aRes1, aSolCriteria, aTouchCriteria;
746 Standard_Real aT1A, aT1B, aT1Min, aT2Min;
747 GeomAPI_ProjectPointOnCurve aProjPC;
748 IntTools_SequenceOfRanges aRanges;
750 aDMin = Precision::Infinite();
751 aSolCriteria = 5.e-16;
752 aTouchCriteria = 5.e-13;
753 Standard_Boolean bTouch = Standard_False;
754 Standard_Boolean bTouchConfirm = Standard_False;
756 aRes1 = Resolution(myCurve1.Curve().Curve(),
757 myCurve1.GetType(), myResCoeff1, myTol);
759 aNbS = SplitRangeOnSegments(aT11, aT12, 3*aRes1, aNbS, aRanges);
761 aProjPC.Init(myGeom2, aT21, aT22);
763 Standard_Real aT11Touch = aT11, aT12Touch = aT12;
764 Standard_Real aT21Touch = aT21, aT22Touch = aT22;
765 Standard_Boolean isSolFound = Standard_False;
766 for (i = 1; i <= aNbS; ++i) {
767 const IntTools_Range& aR1 = aRanges(i);
768 aR1.Range(aT1A, aT1B);
771 iErr = FindDistPC(aT1A, aT1B, myGeom1, aSolCriteria, myPTol1,
772 aProjPC, aD, aT1Min, aT2Min, Standard_False);
778 isSolFound = Standard_True;
781 if (aD < aTouchCriteria) {
785 bTouchConfirm = Standard_True;
790 bTouch = Standard_True;
795 if (!isSolFound || bTouchConfirm)
797 aT1 = (aT11Touch + aT12Touch) * 0.5;
798 iErr = DistPC(aT1, myGeom1, aSolCriteria, aProjPC, aD, aT2, -1);
800 aT2 = (aT21Touch + aT22Touch) * 0.5;
805 //=======================================================================
806 //function : ComputeLineLine
808 //=======================================================================
809 void IntTools_EdgeEdge::ComputeLineLine()
811 Standard_Boolean IsParallel, IsCoincide;
812 Standard_Real aSin, aCos, aAng, aTol;
813 Standard_Real aT1, aT2, aT11, aT12, aT21, aT22;
817 IntTools_CommonPrt aCommonPrt;
819 IsParallel = Standard_False;
820 IsCoincide = Standard_False;
822 aL1 = myCurve1.Line();
823 aL2 = myCurve2.Line();
824 aD1 = aL1.Position().Direction();
825 aD2 = aL2.Position().Direction();
826 myRange1.Range(aT11, aT12);
827 myRange2.Range(aT21, aT22);
829 aCommonPrt.SetEdge1(myEdge1);
830 aCommonPrt.SetEdge2(myEdge2);
833 aAng = (aCos >= 0.) ? 2.*(1. - aCos) : 2.*(1. + aCos);
835 if(aAng <= Precision::Angular()) {
836 IsParallel = Standard_True;
837 if(aL1.SquareDistance(aL2.Location()) <= aTol) {
838 IsCoincide = Standard_True;
839 aP11 = ElCLib::Value(aT11, aL1);
840 aP12 = ElCLib::Value(aT12, aL1);
844 aP11 = ElCLib::Value(aT11, aL1);
845 aP12 = ElCLib::Value(aT12, aL1);
846 if(aL2.SquareDistance(aP11) <= aTol && aL2.SquareDistance(aP12) <= aTol) {
847 IsCoincide = Standard_True;
852 Standard_Real t21, t22;
854 t21 = ElCLib::Parameter(aL2, aP11);
855 t22 = ElCLib::Parameter(aL2, aP12);
856 if((t21 > aT22 && t22 > aT22) || (t21 < aT21 && t22 < aT21)) {
869 aCommonPrt.SetRange1(aT11, aT12);
870 aCommonPrt.SetAllNullFlag(Standard_True);
871 aCommonPrt.AppendRange2(t21, t22);
874 aCommonPrt.SetRange1(aT11, aT12 - (t22 - aT22));
875 aCommonPrt.AppendRange2(t21, aT22);
879 aCommonPrt.SetRange1(aT11 + (aT21 - t21), aT12);
880 aCommonPrt.AppendRange2(aT21, t22);
882 aCommonPrt.SetType(TopAbs_EDGE);
883 myCommonParts.Append(aCommonPrt);
892 TopoDS_Iterator aIt1, aIt2;
893 aIt1.Initialize(myEdge1);
894 for (; aIt1.More(); aIt1.Next()) {
895 const TopoDS_Shape& aV1 = aIt1.Value();
896 aIt2.Initialize(myEdge2);
897 for (; aIt2.More(); aIt2.Next()) {
898 const TopoDS_Shape& aV2 = aIt2.Value();
899 if (aV2.IsSame(aV1)) {
906 aSin = 1. - aCos*aCos;
907 gp_Pnt O1 = aL1.Location();
908 gp_Pnt O2 = aL2.Location();
909 gp_Vec O1O2 (O1, O2);
911 aT2 = (aD1.XYZ()*(O1O2.Dot(aD1))-(O1O2.XYZ())).Dot(aD2.XYZ());
914 if(aT2 < aT21 || aT2 > aT22) {
918 gp_Pnt aP2(ElCLib::Value(aT2, aL2));
919 aT1 = (gp_Vec(O1, aP2)).Dot(aD1);
921 if(aT1 < aT11 || aT1 > aT12) {
925 gp_Pnt aP1(ElCLib::Value(aT1, aL1));
926 Standard_Real aDist = aP1.SquareDistance(aP2);
932 // compute correct range on the edges
933 Standard_Real anAngle, aDt1, aDt2;
935 anAngle = aD1.Angle(aD2);
937 aDt1 = IntTools_Tools::ComputeIntRange(myTol1, myTol2, anAngle);
938 aDt2 = IntTools_Tools::ComputeIntRange(myTol2, myTol1, anAngle);
940 aCommonPrt.SetRange1(aT1 - aDt1, aT1 + aDt1);
941 aCommonPrt.AppendRange2(aT2 - aDt2, aT2 + aDt2);
942 aCommonPrt.SetType(TopAbs_VERTEX);
943 aCommonPrt.SetVertexParameter1(aT1);
944 aCommonPrt.SetVertexParameter2(aT2);
945 myCommonParts.Append(aCommonPrt);
948 //=======================================================================
949 //function : IsIntersection
951 //=======================================================================
952 Standard_Boolean IntTools_EdgeEdge::IsIntersection(const Standard_Real aT11,
953 const Standard_Real aT12,
954 const Standard_Real aT21,
955 const Standard_Real aT22)
957 Standard_Boolean bRet;
958 gp_Pnt aP11, aP12, aP21, aP22;
959 gp_Vec aV11, aV12, aV21, aV22;
960 Standard_Real aD11_21, aD11_22, aD12_21, aD12_22, aCriteria, aCoef;
961 Standard_Boolean bSmall_11_21, bSmall_11_22, bSmall_12_21, bSmall_12_22;
963 bRet = Standard_True;
965 if (((aT12 - aT11) > aCoef*myRes1) && ((aT22 - aT21) > aCoef*myRes2)) {
968 Standard_Real aTRMin = Min((aT12 - aT11)/myRes1, (aT22 - aT21)/myRes2);
969 aCoef = aTRMin / 100.;
974 aCriteria = aCoef * myTol;
975 aCriteria *= aCriteria;
977 myGeom1->D1(aT11, aP11, aV11);
978 myGeom1->D1(aT12, aP12, aV12);
979 myGeom2->D1(aT21, aP21, aV21);
980 myGeom2->D1(aT22, aP22, aV22);
982 aD11_21 = aP11.SquareDistance(aP21);
983 aD11_22 = aP11.SquareDistance(aP22);
984 aD12_21 = aP12.SquareDistance(aP21);
985 aD12_22 = aP12.SquareDistance(aP22);
987 bSmall_11_21 = aD11_21 < aCriteria;
988 bSmall_11_22 = aD11_22 < aCriteria;
989 bSmall_12_21 = aD12_21 < aCriteria;
990 bSmall_12_22 = aD12_22 < aCriteria;
992 if ((bSmall_11_21 && bSmall_12_22) ||
993 (bSmall_11_22 && bSmall_12_21)) {
998 Standard_Real anAngleCriteria;
999 Standard_Real anAngle1, anAngle2;
1001 anAngleCriteria = 5.e-3;
1002 if (bSmall_11_21 && bSmall_12_22) {
1003 anAngle1 = aV11.Angle(aV21);
1004 anAngle2 = aV12.Angle(aV22);
1006 anAngle1 = aV11.Angle(aV22);
1007 anAngle2 = aV12.Angle(aV21);
1010 if (((anAngle1 < anAngleCriteria) || ((M_PI - anAngle1) < anAngleCriteria)) ||
1011 ((anAngle2 < anAngleCriteria) || ((M_PI - anAngle2) < anAngleCriteria))) {
1012 GeomAPI_ProjectPointOnCurve aProjPC;
1013 Standard_Integer iErr;
1014 Standard_Real aD, aT1Min, aT2Min;
1016 aD = Precision::Infinite();
1017 aProjPC.Init(myGeom2, aT21, aT22);
1018 iErr = FindDistPC(aT11, aT12, myGeom1, myTol, myRes1,
1019 aProjPC, aD, aT1Min, aT2Min, Standard_False);
1026 //=======================================================================
1027 //function : CheckCoincidence
1029 //=======================================================================
1030 Standard_Integer IntTools_EdgeEdge::CheckCoincidence(const Standard_Real aT11,
1031 const Standard_Real aT12,
1032 const Standard_Real aT21,
1033 const Standard_Real aT22,
1034 const Standard_Real theCriteria,
1035 const Standard_Real theCurveRes1)
1037 Standard_Integer iErr, aNb, aNb1, i;
1038 Standard_Real aT1A, aT1B, aT1max, aT2max, aDmax;
1039 GeomAPI_ProjectPointOnCurve aProjPC;
1040 IntTools_SequenceOfRanges aRanges;
1044 aProjPC.Init(myGeom2, aT21, aT22);
1046 // 1. Express evaluation
1047 aNb = 10; // Number of intervals on the curve #1
1048 aNb1 = SplitRangeOnSegments(aT11, aT12, theCurveRes1, aNb, aRanges);
1049 for (i = 1; i < aNb1; ++i) {
1050 const IntTools_Range& aR1 = aRanges(i);
1051 aR1.Range(aT1A, aT1B);
1053 iErr = DistPC(aT1B, myGeom1, theCriteria, aProjPC, aDmax, aT2max);
1059 // if the ranges in aRanges are less than theCurveRes1,
1060 // there is no need to do step 2 (deep evaluation)
1065 // 2. Deep evaluation
1066 for (i = 2; i < aNb1; ++i) {
1067 const IntTools_Range& aR1 = aRanges(i);
1068 aR1.Range(aT1A, aT1B);
1070 iErr = FindDistPC(aT1A, aT1B, myGeom1, theCriteria, theCurveRes1,
1071 aProjPC, aDmax, aT1max, aT2max);
1077 // iErr == 0 - the patches are coincided
1078 // iErr == 1 - a point from aC1 can not be projected on aC2
1079 // iErr == 2 - the distance is too big
1083 //=======================================================================
1084 //function : FindDistPC
1086 //=======================================================================
1087 Standard_Integer FindDistPC(const Standard_Real aT1A,
1088 const Standard_Real aT1B,
1089 const Handle(Geom_Curve)& theC1,
1090 const Standard_Real theCriteria,
1091 const Standard_Real theEps,
1092 GeomAPI_ProjectPointOnCurve& theProjPC,
1093 Standard_Real& aDmax,
1094 Standard_Real& aT1max,
1095 Standard_Real& aT2max,
1096 const Standard_Boolean bMaxDist)
1098 Standard_Integer iErr, iC;
1099 Standard_Real aGS, aXP, aA, aB, aXL, aYP, aYL, aT2P, aT2L;
1101 iC = bMaxDist ? 1 : -1;
1103 aT1max = aT2max = 0.; // silence GCC warning
1105 aGS = 0.6180339887498948482045868343656;// =0.5*(1.+sqrt(5.))-1.;
1110 iErr = DistPC(aA, theC1, theCriteria, theProjPC,
1111 aYP, aT2P, aDmax, aT1max, aT2max, iC);
1116 iErr = DistPC(aB, theC1, theCriteria, theProjPC,
1117 aYL, aT2L, aDmax, aT1max, aT2max, iC);
1122 aXP = aA + (aB - aA)*aGS;
1123 aXL = aB - (aB - aA)*aGS;
1125 iErr = DistPC(aXP, theC1, theCriteria, theProjPC,
1126 aYP, aT2P, aDmax, aT1max, aT2max, iC);
1131 iErr = DistPC(aXL, theC1, theCriteria, theProjPC,
1132 aYL, aT2L, aDmax, aT1max, aT2max, iC);
1137 Standard_Real anEps = Max(theEps, Epsilon(Max(Abs(aA), Abs(aB))) * 10.);
1139 if (iC*(aYP - aYL) > 0) {
1143 aXP = aA + (aB - aA)*aGS;
1144 iErr = DistPC(aXP, theC1, theCriteria, theProjPC,
1145 aYP, aT2P, aDmax, aT1max, aT2max, iC);
1151 aXL = aB - (aB - aA)*aGS;
1152 iErr = DistPC(aXL, theC1, theCriteria, theProjPC,
1153 aYL, aT2L, aDmax, aT1max, aT2max, iC);
1157 if ((iErr == 2) && !bMaxDist) {
1158 aXP = (aA + aB) * 0.5;
1159 DistPC(aXP, theC1, theCriteria, theProjPC,
1160 aYP, aT2P, aDmax, aT1max, aT2max, iC);
1165 if ((aB - aA) < anEps) {
1172 //=======================================================================
1175 //=======================================================================
1176 Standard_Integer DistPC(const Standard_Real aT1,
1177 const Handle(Geom_Curve)& theC1,
1178 const Standard_Real theCriteria,
1179 GeomAPI_ProjectPointOnCurve& theProjPC,
1182 Standard_Real& aDmax,
1183 Standard_Real& aT1max,
1184 Standard_Real& aT2max,
1185 const Standard_Integer iC)
1187 Standard_Integer iErr;
1189 iErr = DistPC(aT1, theC1, theCriteria, theProjPC, aD, aT2, iC);
1194 if (iC*(aD - aDmax) > 0) {
1202 //=======================================================================
1205 //=======================================================================
1206 Standard_Integer DistPC(const Standard_Real aT1,
1207 const Handle(Geom_Curve)& theC1,
1208 const Standard_Real theCriteria,
1209 GeomAPI_ProjectPointOnCurve& theProjPC,
1212 const Standard_Integer iC)
1214 Standard_Integer iErr, aNbP2;
1218 theC1->D0(aT1, aP1);
1220 theProjPC.Perform(aP1);
1221 aNbP2 = theProjPC.NbPoints();
1223 iErr = 1;// the point from aC1 can not be projected on aC2
1227 aD = theProjPC.LowerDistance();
1228 aT2 = theProjPC.LowerDistanceParameter();
1229 if (iC*(aD - theCriteria) > 0) {
1230 iErr = 2;// the distance is too big or small
1236 //=======================================================================
1237 //function : SplitRangeOnSegments
1239 //=======================================================================
1240 Standard_Integer SplitRangeOnSegments(const Standard_Real aT1,
1241 const Standard_Real aT2,
1242 const Standard_Real theResolution,
1243 const Standard_Integer theNbSeg,
1244 IntTools_SequenceOfRanges& theSegments)
1246 Standard_Real aDiff = aT2 - aT1;
1247 if (aDiff < theResolution || theNbSeg == 1) {
1248 theSegments.Append(IntTools_Range(aT1, aT2));
1252 Standard_Real aDt, aT1x, aT2x, aSeg;
1253 Standard_Integer aNbSegments, i;
1255 aNbSegments = theNbSeg;
1256 aDt = aDiff / aNbSegments;
1257 if (aDt < theResolution) {
1258 aSeg = aDiff / theResolution;
1259 aNbSegments = Standard_Integer(aSeg) + 1;
1260 aDt = aDiff / aNbSegments;
1264 for (i = 1; i < aNbSegments; ++i) {
1267 IntTools_Range aR(aT1x, aT2x);
1268 theSegments.Append(aR);
1273 IntTools_Range aR(aT1x, aT2);
1274 theSegments.Append(aR);
1279 //=======================================================================
1280 //function : BndBuildBox
1282 //=======================================================================
1283 void BndBuildBox(const BRepAdaptor_Curve& theBAC,
1284 const Standard_Real aT1,
1285 const Standard_Real aT2,
1286 const Standard_Real theTol,
1290 BndLib_Add3dCurve::Add(theBAC, aT1, aT2, theTol, aB);
1294 //=======================================================================
1295 //function : PointBoxDistance
1297 //=======================================================================
1298 Standard_Real PointBoxDistance(const Bnd_Box& aB,
1301 Standard_Real aPCoord[3];
1302 Standard_Real aBMinCoord[3], aBMaxCoord[3];
1303 Standard_Real aDist, aR1, aR2;
1306 aP.Coord(aPCoord[0], aPCoord[1], aPCoord[2]);
1307 aB.Get(aBMinCoord[0], aBMinCoord[1], aBMinCoord[2],
1308 aBMaxCoord[0], aBMaxCoord[1], aBMaxCoord[2]);
1311 for (i = 0; i < 3; ++i) {
1312 aR1 = aBMinCoord[i] - aPCoord[i];
1318 aR2 = aPCoord[i] - aBMaxCoord[i];
1324 aDist = Sqrt(aDist);
1328 //=======================================================================
1329 //function : TypeToInteger
1331 //=======================================================================
1332 Standard_Integer TypeToInteger(const GeomAbs_CurveType theCType)
1334 Standard_Integer iRet;
1340 case GeomAbs_Hyperbola:
1341 case GeomAbs_Parabola:
1344 case GeomAbs_Circle:
1345 case GeomAbs_Ellipse:
1348 case GeomAbs_BezierCurve:
1349 case GeomAbs_BSplineCurve:
1359 //=======================================================================
1360 //function : ResolutionCoeff
1362 //=======================================================================
1363 Standard_Real ResolutionCoeff(const BRepAdaptor_Curve& theBAC,
1364 const IntTools_Range& theRange)
1366 Standard_Real aResCoeff = 0.;
1368 const Handle(Geom_Curve)& aCurve = theBAC.Curve().Curve();
1369 const GeomAbs_CurveType aCurveType = theBAC.GetType();
1371 switch (aCurveType) {
1372 case GeomAbs_Circle :
1373 aResCoeff = 1. / (2 * Handle(Geom_Circle)::DownCast (aCurve)->Circ().Radius());
1375 case GeomAbs_Ellipse :
1376 aResCoeff = 1. / Handle(Geom_Ellipse)::DownCast (aCurve)->MajorRadius();
1378 case GeomAbs_OffsetCurve : {
1379 const Handle(Geom_OffsetCurve)& anOffsetCurve = Handle(Geom_OffsetCurve)::DownCast(aCurve);
1380 const Handle(Geom_Curve)& aBasisCurve = anOffsetCurve->BasisCurve();
1381 GeomAdaptor_Curve aGBasisCurve(aBasisCurve);
1382 const GeomAbs_CurveType aBCType = aGBasisCurve.GetType();
1383 if (aBCType == GeomAbs_Line) {
1386 else if (aBCType == GeomAbs_Circle) {
1387 aResCoeff = 1. / (2 * (anOffsetCurve->Offset() + aGBasisCurve.Circle().Radius()));
1390 else if (aBCType == GeomAbs_Ellipse) {
1391 aResCoeff = 1. / (anOffsetCurve->Offset() + aGBasisCurve.Ellipse().MajorRadius());
1395 Standard_FALLTHROUGH
1396 case GeomAbs_Hyperbola :
1397 case GeomAbs_Parabola :
1398 case GeomAbs_OtherCurve :{
1399 Standard_Real k, kMin, aDist, aDt, aT1, aT2, aT;
1400 Standard_Integer aNbP, i;
1404 theRange.Range(aT1, aT2);
1405 aDt = (aT2 - aT1) / aNbP;
1409 theBAC.D0(aT1, aP1);
1410 for (i = 1; i <= aNbP; ++i) {
1413 aDist = aP1.Distance(aP2);
1431 //=======================================================================
1432 //function : Resolution
1434 //=======================================================================
1435 Standard_Real Resolution(const Handle(Geom_Curve)& theCurve,
1436 const GeomAbs_CurveType theCurveType,
1437 const Standard_Real theResCoeff,
1438 const Standard_Real theR3D)
1442 switch (theCurveType) {
1446 case GeomAbs_Circle: {
1447 Standard_Real aDt = theResCoeff * theR3D;
1448 aRes = (aDt <= 1.) ? 2*ASin(aDt) : 2*M_PI;
1451 case GeomAbs_BezierCurve:
1452 Handle(Geom_BezierCurve)::DownCast (theCurve)->Resolution(theR3D, aRes);
1454 case GeomAbs_BSplineCurve:
1455 Handle(Geom_BSplineCurve)::DownCast (theCurve)->Resolution(theR3D, aRes);
1457 case GeomAbs_OffsetCurve: {
1458 const Handle(Geom_Curve)& aBasisCurve =
1459 Handle(Geom_OffsetCurve)::DownCast(theCurve)->BasisCurve();
1460 const GeomAbs_CurveType aBCType = GeomAdaptor_Curve(aBasisCurve).GetType();
1461 if (aBCType == GeomAbs_Line) {
1465 else if (aBCType == GeomAbs_Circle) {
1466 Standard_Real aDt = theResCoeff * theR3D;
1467 aRes = (aDt <= 1.) ? 2*ASin(aDt) : 2*M_PI;
1471 Standard_FALLTHROUGH
1473 aRes = theResCoeff * theR3D;
1480 //=======================================================================
1481 //function : CurveDeflection
1483 //=======================================================================
1484 Standard_Real CurveDeflection(const BRepAdaptor_Curve& theBAC,
1485 const IntTools_Range& theRange)
1487 Standard_Real aDt, aT, aT1, aT2, aDefl;
1488 Standard_Integer i, aNbP;
1494 theRange.Range(aT1, aT2);
1495 aDt = (aT2 - aT1) / aNbP;
1498 theBAC.D1(aT1, aP, aV1);
1499 for (i = 1; i <= aNbP; ++i) {
1501 theBAC.D1(aT, aP, aV2);
1502 if (aV1.Magnitude() > gp::Resolution() &&
1503 aV2.Magnitude() > gp::Resolution()) {
1504 gp_Dir aD1(aV1), aD2(aV2);
1505 aDefl += aD1.Angle(aD2);
1513 //=======================================================================
1514 //function : IsClosed
1516 //=======================================================================
1517 Standard_Boolean IsClosed(const Handle(Geom_Curve)& theCurve,
1518 const Standard_Real aT1,
1519 const Standard_Real aT2,
1520 const Standard_Real theTol,
1521 const Standard_Real theRes)
1523 if (Abs(aT1 - aT2) < theRes)
1525 return Standard_False;
1529 theCurve->D0(aT1, aP1);
1530 theCurve->D0(aT2, aP2);
1532 Standard_Real aD = aP1.Distance(aP2);