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_Circle.hxx>
25 #include <Geom_Curve.hxx>
26 #include <Geom_Ellipse.hxx>
27 #include <GeomAPI_ProjectPointOnCurve.hxx>
30 #include <IntTools_CommonPrt.hxx>
31 #include <IntTools_EdgeEdge.hxx>
32 #include <IntTools_Range.hxx>
33 #include <TopoDS_Edge.hxx>
34 #include <TopoDS_Iterator.hxx>
37 void BndBuildBox(const BRepAdaptor_Curve& theBAC,
38 const Standard_Real aT1,
39 const Standard_Real aT2,
40 const Standard_Real theTol,
43 Standard_Real PointBoxDistance(const Bnd_Box& aB,
46 Standard_Integer SplitRangeOnSegments(const Standard_Real aT1,
47 const Standard_Real aT2,
48 const Standard_Real theResolution,
49 const Standard_Integer theNbSeg,
50 IntTools_SequenceOfRanges& theSegments);
52 Standard_Integer DistPC(const Standard_Real aT1,
53 const Handle(Geom_Curve)& theC1,
54 const Standard_Real theCriteria,
55 GeomAPI_ProjectPointOnCurve& theProjector,
58 const Standard_Integer iC = 1);
60 Standard_Integer DistPC(const Standard_Real aT1,
61 const Handle(Geom_Curve)& theC1,
62 const Standard_Real theCriteria,
63 GeomAPI_ProjectPointOnCurve& theProjector,
67 Standard_Real& aT1max,
68 Standard_Real& aT2max,
69 const Standard_Integer iC = 1);
71 Standard_Integer FindDistPC(const Standard_Real aT1A,
72 const Standard_Real aT1B,
73 const Handle(Geom_Curve)& theC1,
74 const Standard_Real theCriteria,
75 const Standard_Real theEps,
76 GeomAPI_ProjectPointOnCurve& theProjector,
78 Standard_Real& aT1max,
79 Standard_Real& aT2max,
80 const Standard_Boolean bMaxDist = Standard_True);
82 Standard_Real ResolutionCoeff(const BRepAdaptor_Curve& theBAC,
83 const IntTools_Range& theRange);
85 Standard_Real Resolution(const Handle(Geom_Curve)& theCurve,
86 const GeomAbs_CurveType theCurveType,
87 const Standard_Real theResCoeff,
88 const Standard_Real theR3D);
90 Standard_Real CurveDeflection(const BRepAdaptor_Curve& theBAC,
91 const IntTools_Range& theRange);
93 Standard_Integer IsClosed(const Handle(Geom_Curve)& theCurve,
94 const Standard_Real aT1,
95 const Standard_Real aT2,
96 const Standard_Real theTol,
97 const Standard_Real theRes);
99 Standard_Integer TypeToInteger(const GeomAbs_CurveType theCType);
101 //=======================================================================
104 //=======================================================================
105 void IntTools_EdgeEdge::Prepare()
107 GeomAbs_CurveType aCT1, aCT2;
108 Standard_Integer iCT1, iCT2;
110 myCurve1.Initialize(myEdge1);
111 myCurve2.Initialize(myEdge2);
113 if (myRange1.First() == 0. && myRange1.Last() == 0.) {
114 myRange1.SetFirst(myCurve1.FirstParameter());
115 myRange1.SetLast (myCurve1.LastParameter());
118 if (myRange2.First() == 0. && myRange2.Last() == 0.) {
119 myRange2.SetFirst(myCurve2.FirstParameter());
120 myRange2.SetLast (myCurve2.LastParameter());
123 aCT1 = myCurve1.GetType();
124 aCT2 = myCurve2.GetType();
126 iCT1 = TypeToInteger(aCT1);
127 iCT2 = TypeToInteger(aCT2);
132 Standard_Real aC1, aC2;
134 aC2 = CurveDeflection(myCurve2, myRange2);
135 aC1 = (aC2 > Precision::Confusion()) ?
136 CurveDeflection(myCurve1, myRange1) : 1.;
145 TopoDS_Edge tmpE = myEdge1;
149 BRepAdaptor_Curve tmpC = myCurve1;
153 IntTools_Range tmpR = myRange1;
157 mySwap = Standard_True;
160 myTol1 = myCurve1.Tolerance();
161 myTol2 = myCurve2.Tolerance();
162 myTol = myTol1 + myTol2;
164 if (iCT1 != 0 || iCT2 != 0) {
165 Standard_Real f, l, aTM;
167 myGeom1 = BRep_Tool::Curve(myEdge1, f, l);
168 myGeom2 = BRep_Tool::Curve(myEdge2, f, l);
170 myResCoeff1 = ResolutionCoeff(myCurve1, myRange1);
171 myResCoeff2 = ResolutionCoeff(myCurve2, myRange2);
173 myRes1 = Resolution(myCurve1.Curve().Curve(), myCurve1.GetType(), myResCoeff1, myTol1);
174 myRes2 = Resolution(myCurve2.Curve().Curve(), myCurve2.GetType(), myResCoeff2, myTol2);
177 aTM = Max(fabs(myRange1.First()), fabs(myRange1.Last()));
179 myPTol1 = 5.e-16 * aTM;
183 aTM = Max(fabs(myRange2.First()), fabs(myRange2.Last()));
185 myPTol2 = 5.e-16 * aTM;
190 //=======================================================================
193 //=======================================================================
194 void IntTools_EdgeEdge::Perform()
205 //3.1. Check Line/Line case
206 if (myCurve1.GetType() == GeomAbs_Line &&
207 myCurve2.GetType() == GeomAbs_Line) {
212 IntTools_SequenceOfRanges aRanges1, aRanges2;
214 //3.2. Find ranges containig solutions
215 Standard_Boolean bSplit2;
216 FindSolutions(aRanges1, aRanges2, bSplit2);
218 //4. Merge solutions and save common parts
219 MergeSolutions(aRanges1, aRanges2, bSplit2);
222 //=======================================================================
223 //function : FindSolutions
225 //=======================================================================
226 void IntTools_EdgeEdge::FindSolutions(IntTools_SequenceOfRanges& theRanges1,
227 IntTools_SequenceOfRanges& theRanges2,
228 Standard_Boolean& bSplit2)
230 Standard_Boolean bIsClosed2;
231 Standard_Real aT11, aT12, aT21, aT22;
234 bSplit2 = Standard_False;
235 myRange1.Range(aT11, aT12);
236 myRange2.Range(aT21, aT22);
238 bIsClosed2 = IsClosed(myGeom2, aT21, aT22, myTol2, myRes2);
242 BndBuildBox(myCurve1, aT11, aT12, myTol1, aB1);
244 gp_Pnt aP = myGeom2->Value(aT21);
245 bIsClosed2 = !aB1.IsOut(aP);
249 BndBuildBox(myCurve2, aT21, aT22, myTol2, aB2);
250 FindSolutions(myRange1, myRange2, aB2, theRanges1, theRanges2);
254 if (!CheckCoincidence(aT11, aT12, aT21, aT22, myTol, myRes1)) {
255 theRanges1.Append(myRange1);
256 theRanges2.Append(myRange2);
260 Standard_Integer i, j, aNb1, aNb2;
261 IntTools_SequenceOfRanges aSegments1, aSegments2;
263 aNb1 = IsClosed(myGeom1, aT11, aT12, myTol1, myRes1) ? 2 : 1;
266 aNb1 = SplitRangeOnSegments(aT11, aT12, myRes1, aNb1, aSegments1);
267 aNb2 = SplitRangeOnSegments(aT21, aT22, myRes2, aNb2, aSegments2);
269 for (i = 1; i <= aNb1; ++i) {
270 const IntTools_Range& aR1 = aSegments1(i);
271 for (j = 1; j <= aNb2; ++j) {
272 const IntTools_Range& aR2 = aSegments2(j);
273 BndBuildBox(myCurve2, aR2.First(), aR2.Last(), myTol2, aB2);
274 FindSolutions(aR1, aR2, aB2, theRanges1, theRanges2);
281 //=======================================================================
282 //function : FindSolutions
284 //=======================================================================
285 void IntTools_EdgeEdge::FindSolutions(const IntTools_Range& theR1,
286 const IntTools_Range& theR2,
287 const Bnd_Box& theBox2,
288 IntTools_SequenceOfRanges& theRanges1,
289 IntTools_SequenceOfRanges& theRanges2)
291 Standard_Boolean bOut, bStop, bThin;
292 Standard_Real aT11, aT12, aT21, aT22;
293 Standard_Real aTB11, aTB12, aTB21, aTB22;
294 Standard_Real aSmallStep1, aSmallStep2;
295 Standard_Integer iCom;
298 theR1.Range(aT11, aT12);
299 theR2.Range(aT21, aT22);
303 bThin = Standard_False;
304 bStop = Standard_False;
313 //1. Build box for first edge and find parameters
314 // of the second one in that box
315 BndBuildBox(myCurve1, aT11, aT12, myTol1, aB1);
316 bOut = aB1.IsOut(aB2);
321 bThin = ((aT12 - aT11) < myRes1) ||
322 (aB1.IsXThin(myTol) && aB1.IsYThin(myTol) && aB1.IsZThin(myTol));
324 bOut = !FindParameters(myCurve2, aTB21, aTB22, myRes2, myPTol2,
325 myResCoeff2, aB1, aT21, aT22);
330 //2. Build box for second edge and find parameters
331 // of the first one in that box
332 BndBuildBox(myCurve2, aT21, aT22, myTol2, aB2);
333 bOut = aB1.IsOut(aB2);
338 bThin = ((aT22 - aT21) < myRes2) ||
339 (aB2.IsXThin(myTol) && aB2.IsYThin(myTol) && aB2.IsZThin(myTol));
341 bOut = !FindParameters(myCurve1, aTB11, aTB12, myRes1, myPTol1,
342 myResCoeff1, aB2, aT11, aT12);
348 //3. Check if it makes sense to continue
349 aSmallStep1 = (aTB12 - aTB11) / 250.;
350 aSmallStep2 = (aTB22 - aTB21) / 250.;
352 if (aSmallStep1 < myRes1) {
353 aSmallStep1 = myRes1;
355 if (aSmallStep2 < myRes2) {
356 aSmallStep2 = myRes2;
359 if (((aT11 - aTB11) < aSmallStep1) && ((aTB12 - aT12) < aSmallStep1) &&
360 ((aT21 - aTB21) < aSmallStep2) && ((aTB22 - aT22) < aSmallStep2)) {
361 bStop = Standard_True;
372 //check curves for coincidence on the ranges
373 iCom = CheckCoincidence(aT11, aT12, aT21, aT22, myTol, myRes1);
375 bThin = Standard_True;
381 //check intermediate points
382 Standard_Boolean bSol;
385 GeomAPI_ProjectPointOnCurve aProjPC;
387 aT1 = (aT11 + aT12) * .5;
388 myGeom1->D0(aT1, aP1);
390 aProjPC.Init(myGeom2, aT21, aT22);
391 aProjPC.Perform(aP1);
393 if (aProjPC.NbPoints()) {
394 bSol = aProjPC.LowerDistance() <= myTol;
400 aT2 = (aT21 + aT22) * .5;
401 myGeom2->D0(aT2, aP2);
403 bSol = aP1.IsEqual(aP2, myTol);
411 IntTools_Range aR1(aT11, aT12), aR2(aT21, aT22);
413 theRanges1.Append(aR1);
414 theRanges2.Append(aR2);
418 if (!IsIntersection(aT11, aT12, aT21, aT22)) {
422 //split ranges on segments and repeat
423 Standard_Integer i, aNb1;
424 IntTools_SequenceOfRanges aSegments1;
426 IntTools_Range aR2(aT21, aT22);
427 BndBuildBox(myCurve2, aT21, aT22, myTol2, aB2);
429 aNb1 = SplitRangeOnSegments(aT11, aT12, myRes1, 3, aSegments1);
430 for (i = 1; i <= aNb1; ++i) {
431 const IntTools_Range& aR1 = aSegments1(i);
432 FindSolutions(aR1, aR2, aB2, theRanges1, theRanges2);
436 //=======================================================================
437 //function : FindParameters
439 //=======================================================================
440 Standard_Boolean IntTools_EdgeEdge::FindParameters(const BRepAdaptor_Curve& theBAC,
441 const Standard_Real aT1,
442 const Standard_Real aT2,
443 const Standard_Real theRes,
444 const Standard_Real thePTol,
445 const Standard_Real theResCoeff,
446 const Bnd_Box& theCBox,
450 Standard_Boolean bRet;
451 Standard_Integer aC, i, k;
452 Standard_Real aCf, aDiff, aDt, aT, aTB, aTOut, aTIn;
453 Standard_Real aDist, aDistP, aDistTol, aTol;
457 bRet = Standard_False;
458 aCf = 0.6180339887498948482045868343656;// =0.5*(1.+sqrt(5.))/2.;
460 aTol = theBAC.Tolerance();
462 aDistTol = Precision::PConfusion();
466 const Handle(Geom_Curve)& aCurve = theBAC.Curve().Curve();
467 const GeomAbs_CurveType aCurveType = theBAC.GetType();
469 for (i = 0; i < 2; ++i) {
470 aTB = !i ? aT1 : aT2;
471 aT = !i ? aT2 : aTB1;
473 bRet = Standard_False;
475 //looking for the point on the edge which is in the box;
476 while (aC*(aT-aTB) >= 0) {
478 aDist = PointBoxDistance(theCBox, aP);
480 if (fabs(aDist - aDistP) < aDistTol) {
481 aDt = Resolution(aCurve, aCurveType, theResCoeff, (++k)*aDist);
484 aDt = Resolution(aCurve, aCurveType, theResCoeff, aDist);
488 bRet = Standard_True;
496 //edge is out of the box;
507 //one point IN, one point OUT; looking for the bounding point;
509 aTOut = aTB - aC*aDt;
510 aDiff = aTIn - aTOut;
511 while (fabs(aDiff) > thePTol) {
512 aTB = aTOut + aDiff*aCf;
514 if (aCBx.IsOut(aP)) {
519 aDiff = aTIn - aTOut;
531 //=======================================================================
532 //function : MergeSolutions
534 //=======================================================================
535 void IntTools_EdgeEdge::MergeSolutions(const IntTools_SequenceOfRanges& theRanges1,
536 const IntTools_SequenceOfRanges& theRanges2,
537 const Standard_Boolean bSplit2)
539 Standard_Integer aNbCP = theRanges1.Length();
544 IntTools_Range aRi1, aRi2, aRj1, aRj2;
545 Standard_Boolean bCond;
546 Standard_Integer i, j;
547 TopAbs_ShapeEnum aType;
548 Standard_Real aT11, aT12, aT21, aT22;
549 Standard_Real aTi11, aTi12, aTi21, aTi22;
550 Standard_Real aTj11, aTj12, aTj21, aTj22;
551 Standard_Real aRes1, aRes2, dTR1, dTR2;
552 BOPCol_MapOfInteger aMI;
554 aRes1 = Resolution(myCurve1.Curve().Curve(),
555 myCurve1.GetType(), myResCoeff1, myTol);
556 aRes2 = Resolution(myCurve2.Curve().Curve(),
557 myCurve2.GetType(), myResCoeff2, myTol);
559 myRange1.Range(aT11, aT12);
560 myRange2.Range(aT21, aT22);
563 aType = TopAbs_VERTEX;
565 for (i = 1; i <= aNbCP;) {
566 if (aMI.Contains(i)) {
571 aRi1 = theRanges1(i);
572 aRi2 = theRanges2(i);
574 aRi1.Range(aTi11, aTi12);
575 aRi2.Range(aTi21, aTi22);
579 for (j = i+1; j <= aNbCP; ++j) {
580 if (aMI.Contains(j)) {
584 aRj1 = theRanges1(j);
585 aRj2 = theRanges2(j);
587 aRj1.Range(aTj11, aTj12);
588 aRj2.Range(aTj21, aTj22);
590 bCond = (fabs(aTi12 - aTj11) < dTR1) ||
591 (bSplit2 && (fabs(aTj12 - aTi11) < dTR1));
592 if (bCond && bSplit2) {
593 bCond = (fabs((Max(aTi22, aTj22) - Min(aTi21, aTj21)) -
594 ((aTi22 - aTi21) + (aTj22 - aTj21))) < dTR2);
598 aTi11 = Min(aTi11, aTj11);
599 aTi12 = Max(aTi12, aTj12);
600 aTi21 = Min(aTi21, aTj21);
601 aTi22 = Max(aTi22, aTj22);
610 if (((fabs(aT11 - aTi11) < myRes1) && (fabs(aT12 - aTi12) < myRes1)) ||
611 ((fabs(aT21 - aTi21) < myRes2) && (fabs(aT22 - aTi22) < myRes2))) {
613 myCommonParts.Clear();
616 AddSolution(aTi11, aTi12, aTi21, aTi22, aType);
617 if (aType == TopAbs_EDGE) {
627 //=======================================================================
628 //function : AddSolution
630 //=======================================================================
631 void IntTools_EdgeEdge::AddSolution(const Standard_Real aT11,
632 const Standard_Real aT12,
633 const Standard_Real aT21,
634 const Standard_Real aT22,
635 const TopAbs_ShapeEnum theType)
637 IntTools_CommonPrt aCPart;
639 aCPart.SetType(theType);
641 aCPart.SetEdge1(myEdge1);
642 aCPart.SetEdge2(myEdge2);
643 aCPart.SetRange1(aT11, aT12);
644 aCPart.AppendRange2(aT21, aT22);
646 aCPart.SetEdge1(myEdge2);
647 aCPart.SetEdge2(myEdge1);
648 aCPart.SetRange1(aT21, aT22);
649 aCPart.AppendRange2(aT11, aT12);
652 if (theType == TopAbs_VERTEX) {
653 Standard_Real aT1, aT2;
655 FindBestSolution(aT11, aT12, aT21, aT22, aT1, aT2);
658 aCPart.SetVertexParameter1(aT1);
659 aCPart.SetVertexParameter2(aT2);
661 aCPart.SetVertexParameter1(aT2);
662 aCPart.SetVertexParameter2(aT1);
665 myCommonParts.Append(aCPart);
668 //=======================================================================
669 //function : FindBestSolution
671 //=======================================================================
672 void IntTools_EdgeEdge::FindBestSolution(const Standard_Real aT11,
673 const Standard_Real aT12,
674 const Standard_Real aT21,
675 const Standard_Real aT22,
679 Standard_Integer i, aNbS, iErr;
680 Standard_Real aDMin, aD, aRes1, aSolCriteria, aTouchCriteria;
681 Standard_Real aT1A, aT1B, aT1Min, aT2Min;
682 Standard_Real aT1Im, aT2Im, aT1Touch;
683 GeomAPI_ProjectPointOnCurve aProjPC;
684 IntTools_SequenceOfRanges aRanges;
685 Standard_Boolean bTouch;
687 aDMin = Precision::Infinite();
688 aSolCriteria = 5.e-16;
689 aTouchCriteria = 5.e-13;
690 bTouch = Standard_False;
693 aRes1 = Resolution(myCurve1.Curve().Curve(),
694 myCurve1.GetType(), myResCoeff1, myTol);
696 aNbS = SplitRangeOnSegments(aT11, aT12, 3*aRes1, aNbS, aRanges);
698 aProjPC.Init(myGeom2, aT21, aT22);
700 aT1 = (aT11 + aT12) * 0.5;
701 iErr = DistPC(aT1, myGeom1, aSolCriteria, aProjPC, aD, aT2, -1);
703 aT2 = (aT21 + aT22) * 0.5;
709 for (i = 1; i <= aNbS; ++i) {
710 const IntTools_Range& aR1 = aRanges(i);
711 aR1.Range(aT1A, aT1B);
714 iErr = FindDistPC(aT1A, aT1B, myGeom1, aSolCriteria, myPTol1,
715 aProjPC, aD, aT1Min, aT2Min, Standard_False);
723 if (aD < aTouchCriteria) {
725 aT1A = (aT1Touch + aT1Min) * 0.5;
726 iErr = DistPC(aT1A, myGeom1, aTouchCriteria,
727 aProjPC, aD, aT2Min, -1);
728 if (aD > aTouchCriteria) {
736 bTouch = Standard_True;
743 //=======================================================================
744 //function : ComputeLineLine
746 //=======================================================================
747 void IntTools_EdgeEdge::ComputeLineLine()
749 Standard_Boolean IsParallel, IsCoincide;
750 Standard_Real aSin, aCos, aAng, aTol;
751 Standard_Real aT1, aT2, aT11, aT12, aT21, aT22;
755 IntTools_CommonPrt aCommonPrt;
757 IsParallel = Standard_False;
758 IsCoincide = Standard_False;
760 aL1 = myCurve1.Line();
761 aL2 = myCurve2.Line();
762 aD1 = aL1.Position().Direction();
763 aD2 = aL2.Position().Direction();
764 myRange1.Range(aT11, aT12);
765 myRange2.Range(aT21, aT22);
767 aCommonPrt.SetEdge1(myEdge1);
768 aCommonPrt.SetEdge2(myEdge2);
771 aAng = (aCos >= 0.) ? 2.*(1. - aCos) : 2.*(1. + aCos);
773 if(aAng <= Precision::Angular()) {
774 IsParallel = Standard_True;
775 if(aL1.SquareDistance(aL2.Location()) <= aTol) {
776 IsCoincide = Standard_True;
777 aP11 = ElCLib::Value(aT11, aL1);
778 aP12 = ElCLib::Value(aT12, aL1);
782 aP11 = ElCLib::Value(aT11, aL1);
783 aP12 = ElCLib::Value(aT12, aL1);
784 if(aL2.SquareDistance(aP11) <= aTol && aL2.SquareDistance(aP12) <= aTol) {
785 IsCoincide = Standard_True;
790 Standard_Real t21, t22;
792 t21 = ElCLib::Parameter(aL2, aP11);
793 t22 = ElCLib::Parameter(aL2, aP12);
794 if((t21 > aT22 && t22 > aT22) || (t21 < aT21 && t22 < aT21)) {
807 aCommonPrt.SetRange1(aT11, aT12);
808 aCommonPrt.SetAllNullFlag(Standard_True);
809 aCommonPrt.AppendRange2(t21, t22);
812 aCommonPrt.SetRange1(aT11, aT12 - (t22 - aT22));
813 aCommonPrt.AppendRange2(t21, aT22);
817 aCommonPrt.SetRange1(aT11 + (aT21 - t21), aT12);
818 aCommonPrt.AppendRange2(aT21, t22);
820 aCommonPrt.SetType(TopAbs_EDGE);
821 myCommonParts.Append(aCommonPrt);
830 TopoDS_Iterator aIt1, aIt2;
831 aIt1.Initialize(myEdge1);
832 for (; aIt1.More(); aIt1.Next()) {
833 const TopoDS_Shape& aV1 = aIt1.Value();
834 aIt2.Initialize(myEdge2);
835 for (; aIt2.More(); aIt2.Next()) {
836 const TopoDS_Shape& aV2 = aIt2.Value();
837 if (aV2.IsSame(aV1)) {
844 aSin = 1. - aCos*aCos;
845 gp_Pnt O1 = aL1.Location();
846 gp_Pnt O2 = aL2.Location();
847 gp_Vec O1O2 (O1, O2);
849 aT2 = (aD1.XYZ()*(O1O2.Dot(aD1))-(O1O2.XYZ())).Dot(aD2.XYZ());
852 if(aT2 < aT21 || aT2 > aT22) {
856 gp_Pnt aP2(ElCLib::Value(aT2, aL2));
857 aT1 = (gp_Vec(O1, aP2)).Dot(aD1);
859 if(aT1 < aT11 || aT1 > aT12) {
863 gp_Pnt aP1(ElCLib::Value(aT1, aL1));
864 Standard_Real aDist = aP1.SquareDistance(aP2);
870 aCommonPrt.SetRange1(aT1 - myTol, aT1 + myTol);
871 aCommonPrt.AppendRange2(aT2 - myTol, aT2 + myTol);
872 aCommonPrt.SetType(TopAbs_VERTEX);
873 aCommonPrt.SetVertexParameter1(aT1);
874 aCommonPrt.SetVertexParameter2(aT2);
875 myCommonParts.Append(aCommonPrt);
878 //=======================================================================
879 //function : IsIntersection
881 //=======================================================================
882 Standard_Boolean IntTools_EdgeEdge::IsIntersection(const Standard_Real aT11,
883 const Standard_Real aT12,
884 const Standard_Real aT21,
885 const Standard_Real aT22)
887 Standard_Boolean bRet;
888 gp_Pnt aP11, aP12, aP21, aP22;
889 gp_Vec aV11, aV12, aV21, aV22;
890 Standard_Real aD11_21, aD11_22, aD12_21, aD12_22, aCriteria, aCoef;
891 Standard_Boolean bSmall_11_21, bSmall_11_22, bSmall_12_21, bSmall_12_22;
893 bRet = Standard_True;
895 if (((aT12 - aT11) > aCoef*myRes1) && ((aT22 - aT21) > aCoef*myRes2)) {
898 Standard_Real aTRMin = Min((aT12 - aT11)/myRes1, (aT22 - aT21)/myRes2);
899 aCoef = aTRMin / 100.;
904 aCriteria = aCoef * myTol;
905 aCriteria *= aCriteria;
907 myGeom1->D1(aT11, aP11, aV11);
908 myGeom1->D1(aT12, aP12, aV12);
909 myGeom2->D1(aT21, aP21, aV21);
910 myGeom2->D1(aT22, aP22, aV22);
912 aD11_21 = aP11.SquareDistance(aP21);
913 aD11_22 = aP11.SquareDistance(aP22);
914 aD12_21 = aP12.SquareDistance(aP21);
915 aD12_22 = aP12.SquareDistance(aP22);
917 bSmall_11_21 = aD11_21 < aCriteria;
918 bSmall_11_22 = aD11_22 < aCriteria;
919 bSmall_12_21 = aD12_21 < aCriteria;
920 bSmall_12_22 = aD12_22 < aCriteria;
922 if ((bSmall_11_21 && bSmall_12_22) ||
923 (bSmall_11_22 && bSmall_12_21)) {
928 Standard_Real anAngleCriteria;
929 Standard_Real anAngle1, anAngle2;
931 anAngleCriteria = 5.e-3;
932 if (bSmall_11_21 && bSmall_12_22) {
933 anAngle1 = aV11.Angle(aV21);
934 anAngle2 = aV12.Angle(aV22);
936 anAngle1 = aV11.Angle(aV22);
937 anAngle2 = aV12.Angle(aV21);
940 if (((anAngle1 < anAngleCriteria) || ((M_PI - anAngle1) < anAngleCriteria)) ||
941 ((anAngle2 < anAngleCriteria) || ((M_PI - anAngle2) < anAngleCriteria))) {
942 GeomAPI_ProjectPointOnCurve aProjPC;
943 Standard_Integer iErr;
944 Standard_Real aD, aT1Min, aT2Min;
946 aD = Precision::Infinite();
947 aProjPC.Init(myGeom2, aT21, aT22);
948 iErr = FindDistPC(aT11, aT12, myGeom1, myTol, myRes1,
949 aProjPC, aD, aT1Min, aT2Min, Standard_False);
956 //=======================================================================
957 //function : CheckCoincidence
959 //=======================================================================
960 Standard_Integer IntTools_EdgeEdge::CheckCoincidence(const Standard_Real aT11,
961 const Standard_Real aT12,
962 const Standard_Real aT21,
963 const Standard_Real aT22,
964 const Standard_Real theCriteria,
965 const Standard_Real theCurveRes1)
967 Standard_Integer iErr, aNb, aNb1, i;
968 Standard_Real aT1A, aT1B, aT1max, aT2max, aDmax;
969 GeomAPI_ProjectPointOnCurve aProjPC;
970 IntTools_SequenceOfRanges aRanges;
974 aProjPC.Init(myGeom2, aT21, aT22);
976 // 1. Express evaluation
977 aNb = 10; // Number of intervals on the curve #1
978 aNb1 = SplitRangeOnSegments(aT11, aT12, theCurveRes1, aNb, aRanges);
979 for (i = 1; i < aNb1; ++i) {
980 const IntTools_Range& aR1 = aRanges(i);
981 aR1.Range(aT1A, aT1B);
983 iErr = DistPC(aT1B, myGeom1, theCriteria, aProjPC, aDmax, aT2max);
989 // if the ranges in aRanges are less than theCurveRes1,
990 // there is no need to do step 2 (deep evaluation)
995 // 2. Deep evaluation
996 for (i = 2; i < aNb1; ++i) {
997 const IntTools_Range& aR1 = aRanges(i);
998 aR1.Range(aT1A, aT1B);
1000 iErr = FindDistPC(aT1A, aT1B, myGeom1, theCriteria, theCurveRes1,
1001 aProjPC, aDmax, aT1max, aT2max);
1007 // iErr == 0 - the patches are coincided
1008 // iErr == 1 - a point from aC1 can not be projected on aC2
1009 // iErr == 2 - the distance is too big
1013 //=======================================================================
1014 //function : FindDistPC
1016 //=======================================================================
1017 Standard_Integer FindDistPC(const Standard_Real aT1A,
1018 const Standard_Real aT1B,
1019 const Handle(Geom_Curve)& theC1,
1020 const Standard_Real theCriteria,
1021 const Standard_Real theEps,
1022 GeomAPI_ProjectPointOnCurve& theProjPC,
1023 Standard_Real& aDmax,
1024 Standard_Real& aT1max,
1025 Standard_Real& aT2max,
1026 const Standard_Boolean bMaxDist)
1028 Standard_Integer iErr, iC;
1029 Standard_Real aGS, aXP, aA, aB, aXL, aYP, aYL, aT2P, aT2L;
1031 iC = bMaxDist ? 1 : -1;
1034 aGS = 0.6180339887498948482045868343656;// =0.5*(1.+sqrt(5.))-1.;
1039 iErr = DistPC(aA, theC1, theCriteria, theProjPC,
1040 aYP, aT2P, aDmax, aT1max, aT2max, iC);
1045 iErr = DistPC(aB, theC1, theCriteria, theProjPC,
1046 aYL, aT2L, aDmax, aT1max, aT2max, iC);
1051 aXP = aA + (aB - aA)*aGS;
1052 aXL = aB - (aB - aA)*aGS;
1054 iErr = DistPC(aXP, theC1, theCriteria, theProjPC,
1055 aYP, aT2P, aDmax, aT1max, aT2max, iC);
1060 iErr = DistPC(aXL, theC1, theCriteria, theProjPC,
1061 aYL, aT2L, aDmax, aT1max, aT2max, iC);
1067 if (iC*(aYP - aYL) > 0) {
1071 aXP = aA + (aB - aA)*aGS;
1072 iErr = DistPC(aXP, theC1, theCriteria, theProjPC,
1073 aYP, aT2P, aDmax, aT1max, aT2max, iC);
1079 aXL = aB - (aB - aA)*aGS;
1080 iErr = DistPC(aXL, theC1, theCriteria, theProjPC,
1081 aYL, aT2L, aDmax, aT1max, aT2max, iC);
1085 if ((iErr == 2) && !bMaxDist) {
1086 aXP = (aA + aB) * 0.5;
1087 DistPC(aXP, theC1, theCriteria, theProjPC,
1088 aYP, aT2P, aDmax, aT1max, aT2max, iC);
1093 if ((aB - aA) < theEps) {
1100 //=======================================================================
1103 //=======================================================================
1104 Standard_Integer DistPC(const Standard_Real aT1,
1105 const Handle(Geom_Curve)& theC1,
1106 const Standard_Real theCriteria,
1107 GeomAPI_ProjectPointOnCurve& theProjPC,
1110 Standard_Real& aDmax,
1111 Standard_Real& aT1max,
1112 Standard_Real& aT2max,
1113 const Standard_Integer iC)
1115 Standard_Integer iErr;
1117 iErr = DistPC(aT1, theC1, theCriteria, theProjPC, aD, aT2, iC);
1122 if (iC*(aD - aDmax) > 0) {
1130 //=======================================================================
1133 //=======================================================================
1134 Standard_Integer DistPC(const Standard_Real aT1,
1135 const Handle(Geom_Curve)& theC1,
1136 const Standard_Real theCriteria,
1137 GeomAPI_ProjectPointOnCurve& theProjPC,
1140 const Standard_Integer iC)
1142 Standard_Integer iErr, aNbP2;
1146 theC1->D0(aT1, aP1);
1148 theProjPC.Perform(aP1);
1149 aNbP2 = theProjPC.NbPoints();
1151 iErr = 1;// the point from aC1 can not be projected on aC2
1155 aD = theProjPC.LowerDistance();
1156 aT2 = theProjPC.LowerDistanceParameter();
1157 if (iC*(aD - theCriteria) > 0) {
1158 iErr = 2;// the distance is too big or small
1164 //=======================================================================
1165 //function : SplitRangeOnSegments
1167 //=======================================================================
1168 Standard_Integer SplitRangeOnSegments(const Standard_Real aT1,
1169 const Standard_Real aT2,
1170 const Standard_Real theResolution,
1171 const Standard_Integer theNbSeg,
1172 IntTools_SequenceOfRanges& theSegments)
1174 Standard_Real aDiff = aT2 - aT1;
1175 if (aDiff < theResolution || theNbSeg == 1) {
1176 theSegments.Append(IntTools_Range(aT1, aT2));
1180 Standard_Real aDt, aT1x, aT2x, aSeg;
1181 Standard_Integer aNbSegments, i;
1183 aNbSegments = theNbSeg;
1184 aDt = aDiff / aNbSegments;
1185 if (aDt < theResolution) {
1186 aSeg = aDiff / theResolution;
1187 aNbSegments = Standard_Integer(aSeg) + 1;
1188 aDt = aDiff / aNbSegments;
1192 for (i = 1; i < aNbSegments; ++i) {
1195 IntTools_Range aR(aT1x, aT2x);
1196 theSegments.Append(aR);
1201 IntTools_Range aR(aT1x, aT2);
1202 theSegments.Append(aR);
1207 //=======================================================================
1208 //function : BndBuildBox
1210 //=======================================================================
1211 void BndBuildBox(const BRepAdaptor_Curve& theBAC,
1212 const Standard_Real aT1,
1213 const Standard_Real aT2,
1214 const Standard_Real theTol,
1218 BndLib_Add3dCurve::Add(theBAC, aT1, aT2, theTol, aB);
1222 //=======================================================================
1223 //function : PointBoxDistance
1225 //=======================================================================
1226 Standard_Real PointBoxDistance(const Bnd_Box& aB,
1229 Standard_Real aPCoord[3];
1230 Standard_Real aBMinCoord[3], aBMaxCoord[3];
1231 Standard_Real aDist, aR1, aR2;
1234 aP.Coord(aPCoord[0], aPCoord[1], aPCoord[2]);
1235 aB.Get(aBMinCoord[0], aBMinCoord[1], aBMinCoord[2],
1236 aBMaxCoord[0], aBMaxCoord[1], aBMaxCoord[2]);
1239 for (i = 0; i < 3; ++i) {
1240 aR1 = aBMinCoord[i] - aPCoord[i];
1246 aR2 = aPCoord[i] - aBMaxCoord[i];
1252 aDist = Sqrt(aDist);
1256 //=======================================================================
1257 //function : TypeToInteger
1259 //=======================================================================
1260 Standard_Integer TypeToInteger(const GeomAbs_CurveType theCType)
1262 Standard_Integer iRet;
1268 case GeomAbs_Hyperbola:
1269 case GeomAbs_Parabola:
1272 case GeomAbs_Circle:
1273 case GeomAbs_Ellipse:
1276 case GeomAbs_BezierCurve:
1277 case GeomAbs_BSplineCurve:
1287 //=======================================================================
1288 //function : ResolutionCoeff
1290 //=======================================================================
1291 Standard_Real ResolutionCoeff(const BRepAdaptor_Curve& theBAC,
1292 const IntTools_Range& theRange)
1294 Standard_Real aResCoeff;
1296 const Handle(Geom_Curve)& aCurve = theBAC.Curve().Curve();
1297 const GeomAbs_CurveType aCurveType = theBAC.GetType();
1299 switch (aCurveType) {
1300 case GeomAbs_Circle :
1301 aResCoeff = 1. / (2 * Handle(Geom_Circle)::DownCast (aCurve)->Circ().Radius());
1303 case GeomAbs_Ellipse :
1304 aResCoeff = 1. / Handle(Geom_Ellipse)::DownCast (aCurve)->MajorRadius();
1306 case GeomAbs_Hyperbola :
1307 case GeomAbs_Parabola :
1308 case GeomAbs_OtherCurve :{
1309 Standard_Real k, kMin, aDist, aDt, aT1, aT2, aT;
1310 Standard_Integer aNbP, i;
1314 theRange.Range(aT1, aT2);
1315 aDt = (aT2 - aT1) / aNbP;
1319 theBAC.D0(aT1, aP1);
1320 for (i = 1; i <= aNbP; ++i) {
1323 aDist = aP1.Distance(aP2);
1342 //=======================================================================
1343 //function : Resolution
1345 //=======================================================================
1346 Standard_Real Resolution(const Handle(Geom_Curve)& theCurve,
1347 const GeomAbs_CurveType theCurveType,
1348 const Standard_Real theResCoeff,
1349 const Standard_Real theR3D)
1353 switch (theCurveType) {
1357 case GeomAbs_Circle: {
1358 Standard_Real aDt = theResCoeff * theR3D;
1359 aRes = (aDt <= 1.) ? 2*ASin(aDt) : 2*M_PI;
1362 case GeomAbs_BezierCurve:
1363 Handle(Geom_BezierCurve)::DownCast (theCurve)->Resolution(theR3D, aRes);
1365 case GeomAbs_BSplineCurve:
1366 Handle(Geom_BSplineCurve)::DownCast (theCurve)->Resolution(theR3D, aRes);
1369 aRes = theResCoeff * theR3D;
1376 //=======================================================================
1377 //function : CurveDeflection
1379 //=======================================================================
1380 Standard_Real CurveDeflection(const BRepAdaptor_Curve& theBAC,
1381 const IntTools_Range& theRange)
1383 Standard_Real aDt, aT, aT1, aT2, aDefl;
1384 Standard_Integer i, aNbP;
1390 theRange.Range(aT1, aT2);
1391 aDt = (aT2 - aT1) / aNbP;
1394 theBAC.D1(aT1, aP, aV1);
1395 for (i = 1; i <= aNbP; ++i) {
1397 theBAC.D1(aT, aP, aV2);
1398 if (aV1.Magnitude() > gp::Resolution() &&
1399 aV2.Magnitude() > gp::Resolution()) {
1400 gp_Dir aD1(aV1), aD2(aV2);
1401 aDefl += aD1.Angle(aD2);
1409 //=======================================================================
1410 //function : IsClosed
1412 //=======================================================================
1413 Standard_Integer IsClosed(const Handle(Geom_Curve)& theCurve,
1414 const Standard_Real aT1,
1415 const Standard_Real aT2,
1416 const Standard_Real theTol,
1417 const Standard_Real theRes)
1419 Standard_Boolean bClosed;
1423 bClosed = Standard_False;
1424 if (Abs(aT1 - aT2) < theRes) {
1428 theCurve->D0(aT1, aP1);
1429 theCurve->D0(aT2, aP2);
1431 aD = aP1.Distance(aP2);
1432 bClosed = aD < theTol;