1 // Created on: 2000-11-23
2 // Created by: Michael KLOKOV
3 // Copyright (c) 2000-2014 OPEN CASCADE SAS
5 // This file is part of Open CASCADE Technology software library.
7 // This library is free software; you can redistribute it and/or modify it under
8 // the terms of the GNU Lesser General Public License version 2.1 as published
9 // by the Free Software Foundation, with special exception defined in the file
10 // OCCT_LGPL_EXCEPTION.txt. Consult the file LICENSE_LGPL_21.txt included in OCCT
11 // distribution for complete text of the license and disclaimer of any warranty.
13 // Alternatively, this file may be used under the terms of Open CASCADE
14 // commercial license or contractual agreement.
16 #include <IntTools_FaceFace.hxx>
18 #include <Precision.hxx>
20 #include <TColStd_HArray1OfReal.hxx>
21 #include <TColStd_Array1OfReal.hxx>
22 #include <TColStd_Array1OfInteger.hxx>
23 #include <TColStd_SequenceOfReal.hxx>
24 #include <TColStd_ListOfInteger.hxx>
25 #include <TColStd_ListIteratorOfListOfInteger.hxx>
26 #include <TColStd_Array1OfListOfInteger.hxx>
28 #include <gp_Lin2d.hxx>
29 #include <gp_Ax22d.hxx>
30 #include <gp_Circ2d.hxx>
31 #include <gp_Torus.hxx>
32 #include <gp_Cylinder.hxx>
34 #include <Bnd_Box.hxx>
36 #include <TColgp_HArray1OfPnt2d.hxx>
37 #include <TColgp_SequenceOfPnt2d.hxx>
38 #include <TColgp_Array1OfPnt.hxx>
39 #include <TColgp_Array1OfPnt2d.hxx>
41 #include <IntAna_QuadQuadGeo.hxx>
43 #include <IntSurf_PntOn2S.hxx>
44 #include <IntSurf_LineOn2S.hxx>
45 #include <IntSurf_PntOn2S.hxx>
46 #include <IntSurf_ListOfPntOn2S.hxx>
47 #include <IntRes2d_Domain.hxx>
48 #include <ProjLib_Plane.hxx>
50 #include <IntPatch_GLine.hxx>
51 #include <IntPatch_RLine.hxx>
52 #include <IntPatch_WLine.hxx>
53 #include <IntPatch_ALine.hxx>
54 #include <IntPatch_ALineToWLine.hxx>
59 #include <Extrema_ExtCC.hxx>
60 #include <Extrema_POnCurv.hxx>
61 #include <BndLib_AddSurface.hxx>
63 #include <Adaptor3d_SurfacePtr.hxx>
64 #include <Adaptor2d_HLine2d.hxx>
66 #include <GeomAbs_SurfaceType.hxx>
67 #include <GeomAbs_CurveType.hxx>
69 #include <Geom_Surface.hxx>
70 #include <Geom_Line.hxx>
71 #include <Geom_Circle.hxx>
72 #include <Geom_Ellipse.hxx>
73 #include <Geom_Parabola.hxx>
74 #include <Geom_Hyperbola.hxx>
75 #include <Geom_TrimmedCurve.hxx>
76 #include <Geom_BSplineCurve.hxx>
77 #include <Geom_RectangularTrimmedSurface.hxx>
78 #include <Geom_OffsetSurface.hxx>
79 #include <Geom_Curve.hxx>
80 #include <Geom_Conic.hxx>
82 #include <Geom2d_TrimmedCurve.hxx>
83 #include <Geom2d_BSplineCurve.hxx>
84 #include <Geom2d_Line.hxx>
85 #include <Geom2d_Curve.hxx>
86 #include <Geom2d_Circle.hxx>
88 #include <Geom2dAPI_InterCurveCurve.hxx>
89 #include <Geom2dInt_GInter.hxx>
90 #include <Geom2dAdaptor.hxx>
91 #include <GeomAdaptor_Curve.hxx>
92 #include <GeomAdaptor_HSurface.hxx>
93 #include <GeomAdaptor_Surface.hxx>
94 #include <GeomLib_CheckBSplineCurve.hxx>
95 #include <GeomLib_Check2dBSplineCurve.hxx>
97 #include <GeomInt_WLApprox.hxx>
98 #include <GeomProjLib.hxx>
99 #include <GeomAPI_ProjectPointOnSurf.hxx>
100 #include <Geom2dAdaptor_Curve.hxx>
101 #include <TopoDS.hxx>
102 #include <TopoDS_Edge.hxx>
103 #include <TopExp_Explorer.hxx>
105 #include <BRep_Tool.hxx>
106 #include <BRepTools.hxx>
107 #include <BRepAdaptor_Surface.hxx>
109 #include <IntTools_Curve.hxx>
110 #include <IntTools_Tools.hxx>
111 #include <IntTools_Tools.hxx>
112 #include <IntTools_TopolTool.hxx>
113 #include <IntTools_PntOnFace.hxx>
114 #include <IntTools_PntOn2Faces.hxx>
115 #include <IntTools_Context.hxx>
116 #include <IntSurf_ListIteratorOfListOfPntOn2S.hxx>
117 #include <GeomInt.hxx>
119 #include <Approx_CurveOnSurface.hxx>
120 #include <GeomAdaptor.hxx>
121 #include <GeomInt_IntSS.hxx>
124 void RefineVector(gp_Vec2d& aV2D);
125 #ifdef OCCT_DEBUG_DUMPWLINE
127 void DumpWLine(const Handle(IntPatch_WLine)& aWLine);
131 void TolR3d(const TopoDS_Face& ,
135 Handle(Geom_Curve) MakeBSpline (const Handle(IntPatch_WLine)&,
136 const Standard_Integer,
137 const Standard_Integer);
140 void Parameters(const Handle(GeomAdaptor_HSurface)&,
141 const Handle(GeomAdaptor_HSurface)&,
149 void CorrectSurfaceBoundaries(const TopoDS_Face& theFace,
150 const Standard_Real theTolerance,
151 Standard_Real& theumin,
152 Standard_Real& theumax,
153 Standard_Real& thevmin,
154 Standard_Real& thevmax);
157 Standard_Boolean NotUseSurfacesForApprox
158 (const TopoDS_Face& aF1,
159 const TopoDS_Face& aF2,
160 const Handle(IntPatch_WLine)& WL,
161 const Standard_Integer ifprm,
162 const Standard_Integer ilprm);
165 Handle(IntPatch_WLine) ComputePurgedWLine(const Handle(IntPatch_WLine) &theWLine,
166 const Handle(GeomAdaptor_HSurface) &theS1,
167 const Handle(GeomAdaptor_HSurface) &theS2);
170 Handle(Geom2d_BSplineCurve) MakeBSpline2d(const Handle(IntPatch_WLine)& theWLine,
171 const Standard_Integer ideb,
172 const Standard_Integer ifin,
173 const Standard_Boolean onFirst);
176 Standard_Boolean DecompositionOfWLine(const Handle(IntPatch_WLine)& theWLine,
177 const Handle(GeomAdaptor_HSurface)& theSurface1,
178 const Handle(GeomAdaptor_HSurface)& theSurface2,
179 const TopoDS_Face& theFace1,
180 const TopoDS_Face& theFace2,
181 const GeomInt_LineConstructor& theLConstructor,
182 const Standard_Boolean theAvoidLConstructor,
183 IntPatch_SequenceOfLine& theNewLines,
184 Standard_Real& theReachedTol3d,
185 const Handle(IntTools_Context)& );
188 Standard_Boolean ParameterOutOfBoundary(const Standard_Real theParameter,
189 const Handle(Geom_Curve)& theCurve,
190 const TopoDS_Face& theFace1,
191 const TopoDS_Face& theFace2,
192 const Standard_Real theOtherParameter,
193 const Standard_Boolean bIncreasePar,
194 Standard_Real& theNewParameter,
195 const Handle(IntTools_Context)& );
198 Standard_Boolean IsCurveValid(Handle(Geom2d_Curve)& thePCurve);
201 Standard_Boolean IsPointOnBoundary(const Standard_Real theParameter,
202 const Standard_Real theFirstBoundary,
203 const Standard_Real theSecondBoundary,
204 const Standard_Real theResolution,
205 Standard_Boolean& IsOnFirstBoundary);
207 Standard_Boolean FindPoint(const gp_Pnt2d& theFirstPoint,
208 const gp_Pnt2d& theLastPoint,
209 const Standard_Real theUmin,
210 const Standard_Real theUmax,
211 const Standard_Real theVmin,
212 const Standard_Real theVmax,
213 gp_Pnt2d& theNewPoint);
217 Standard_Integer ComputeTangentZones( const Handle(GeomAdaptor_HSurface)& theSurface1,
218 const Handle(GeomAdaptor_HSurface)& theSurface2,
219 const TopoDS_Face& theFace1,
220 const TopoDS_Face& theFace2,
221 Handle(TColgp_HArray1OfPnt2d)& theResultOnS1,
222 Handle(TColgp_HArray1OfPnt2d)& theResultOnS2,
223 Handle(TColStd_HArray1OfReal)& theResultRadius,
224 const Handle(IntTools_Context)& );
227 Standard_Boolean FindPoint(const gp_Pnt2d& theFirstPoint,
228 const gp_Pnt2d& theLastPoint,
229 const Standard_Real theUmin,
230 const Standard_Real theUmax,
231 const Standard_Real theVmin,
232 const Standard_Real theVmax,
233 const gp_Pnt2d& theTanZoneCenter,
234 const Standard_Real theZoneRadius,
235 Handle(GeomAdaptor_HSurface) theGASurface,
236 gp_Pnt2d& theNewPoint);
239 Standard_Boolean IsInsideTanZone(const gp_Pnt2d& thePoint,
240 const gp_Pnt2d& theTanZoneCenter,
241 const Standard_Real theZoneRadius,
242 Handle(GeomAdaptor_HSurface) theGASurface);
245 gp_Pnt2d AdjustByNeighbour(const gp_Pnt2d& theaNeighbourPoint,
246 const gp_Pnt2d& theOriginalPoint,
247 Handle(GeomAdaptor_HSurface) theGASurface);
249 Standard_Boolean ApproxWithPCurves(const gp_Cylinder& theCyl,
250 const gp_Sphere& theSph);
252 static void PerformPlanes(const Handle(GeomAdaptor_HSurface)& theS1,
253 const Handle(GeomAdaptor_HSurface)& theS2,
254 const Standard_Real TolAng,
255 const Standard_Real TolTang,
256 const Standard_Boolean theApprox1,
257 const Standard_Boolean theApprox2,
258 IntTools_SequenceOfCurves& theSeqOfCurve,
259 Standard_Boolean& theTangentFaces);
261 static Standard_Boolean ClassifyLin2d(const Handle(GeomAdaptor_HSurface)& theS,
262 const gp_Lin2d& theLin2d,
263 const Standard_Real theTol,
264 Standard_Real& theP1,
265 Standard_Real& theP2);
268 void ApproxParameters(const Handle(GeomAdaptor_HSurface)& aHS1,
269 const Handle(GeomAdaptor_HSurface)& aHS2,
270 Standard_Integer& iDegMin,
271 Standard_Integer& iNbIter,
272 Standard_Integer& iDegMax);
275 void Tolerances(const Handle(GeomAdaptor_HSurface)& aHS1,
276 const Handle(GeomAdaptor_HSurface)& aHS2,
277 Standard_Real& aTolTang);
280 Standard_Boolean SortTypes(const GeomAbs_SurfaceType aType1,
281 const GeomAbs_SurfaceType aType2);
283 Standard_Integer IndexType(const GeomAbs_SurfaceType aType);
287 Standard_Boolean CheckPCurve(const Handle(Geom2d_Curve)& aPC,
288 const TopoDS_Face& aFace);
291 Standard_Real MaxDistance(const Handle(Geom_Curve)& theC,
292 const Standard_Real aT,
293 GeomAPI_ProjectPointOnSurf& theProjPS);
296 Standard_Real FindMaxDistance(const Handle(Geom_Curve)& theC,
297 const Standard_Real theFirst,
298 const Standard_Real theLast,
299 GeomAPI_ProjectPointOnSurf& theProjPS,
300 const Standard_Real theEps);
303 Standard_Real FindMaxDistance(const Handle(Geom_Curve)& theCurve,
304 const Standard_Real theFirst,
305 const Standard_Real theLast,
306 const TopoDS_Face& theFace,
307 const Handle(IntTools_Context)& theContext);
310 void CorrectPlaneBoundaries(Standard_Real& aUmin,
311 Standard_Real& aUmax,
312 Standard_Real& aVmin,
313 Standard_Real& aVmax);
315 //=======================================================================
318 //=======================================================================
319 IntTools_FaceFace::IntTools_FaceFace()
321 myIsDone=Standard_False;
322 myTangentFaces=Standard_False;
324 myHS1 = new GeomAdaptor_HSurface ();
325 myHS2 = new GeomAdaptor_HSurface ();
328 SetParameters(Standard_True, Standard_True, Standard_True, 1.e-07);
331 //=======================================================================
332 //function : SetContext
334 //=======================================================================
335 void IntTools_FaceFace::SetContext(const Handle(IntTools_Context)& aContext)
339 //=======================================================================
342 //=======================================================================
343 const Handle(IntTools_Context)& IntTools_FaceFace::Context()const
347 //=======================================================================
350 //=======================================================================
351 const TopoDS_Face& IntTools_FaceFace::Face1() const
355 //=======================================================================
358 //=======================================================================
359 const TopoDS_Face& IntTools_FaceFace::Face2() const
363 //=======================================================================
364 //function : TangentFaces
366 //=======================================================================
367 Standard_Boolean IntTools_FaceFace::TangentFaces() const
369 return myTangentFaces;
371 //=======================================================================
374 //=======================================================================
375 const IntTools_SequenceOfPntOn2Faces& IntTools_FaceFace::Points() const
379 //=======================================================================
382 //=======================================================================
383 Standard_Boolean IntTools_FaceFace::IsDone() const
387 //=======================================================================
388 //function : TolReached3d
390 //=======================================================================
391 Standard_Real IntTools_FaceFace::TolReached3d() const
393 return myTolReached3d;
395 //=======================================================================
397 //purpose : return lines of intersection
398 //=======================================================================
399 const IntTools_SequenceOfCurves& IntTools_FaceFace::Lines() const
401 StdFail_NotDone_Raise_if
403 "IntTools_FaceFace::Lines() => myIntersector NOT DONE");
406 //=======================================================================
407 //function : TolReached2d
409 //=======================================================================
410 Standard_Real IntTools_FaceFace::TolReached2d() const
412 return myTolReached2d;
414 // =======================================================================
415 // function: SetParameters
417 // =======================================================================
418 void IntTools_FaceFace::SetParameters(const Standard_Boolean ToApproxC3d,
419 const Standard_Boolean ToApproxC2dOnS1,
420 const Standard_Boolean ToApproxC2dOnS2,
421 const Standard_Real ApproximationTolerance)
423 myApprox = ToApproxC3d;
424 myApprox1 = ToApproxC2dOnS1;
425 myApprox2 = ToApproxC2dOnS2;
426 myTolApprox = ApproximationTolerance;
428 //=======================================================================
431 //=======================================================================
432 void IntTools_FaceFace::SetList(IntSurf_ListOfPntOn2S& aListOfPnts)
434 myListOfPnts = aListOfPnts;
438 static Standard_Boolean isTreatAnalityc(const TopoDS_Face& theF1,
439 const TopoDS_Face& theF2)
441 const Standard_Real Tolang = 1.e-8;
442 const Standard_Real aTolF1=BRep_Tool::Tolerance(theF1);
443 const Standard_Real aTolF2=BRep_Tool::Tolerance(theF2);
444 const Standard_Real aTolSum = aTolF1 + aTolF2;
445 Standard_Real aHigh = 0.0;
447 const BRepAdaptor_Surface aBAS1(theF1), aBAS2(theF2);
448 const GeomAbs_SurfaceType aType1=aBAS1.GetType();
449 const GeomAbs_SurfaceType aType2=aBAS2.GetType();
453 if(aType1 == GeomAbs_Plane)
457 else if(aType2 == GeomAbs_Plane)
463 return Standard_True;
466 if(aType1 == GeomAbs_Cylinder)
468 aS2=aBAS1.Cylinder();
469 const Standard_Real VMin = aBAS1.FirstVParameter();
470 const Standard_Real VMax = aBAS1.LastVParameter();
472 if( Precision::IsNegativeInfinite(VMin) ||
473 Precision::IsPositiveInfinite(VMax))
474 return Standard_True;
478 else if(aType2 == GeomAbs_Cylinder)
480 aS2=aBAS2.Cylinder();
482 const Standard_Real VMin = aBAS2.FirstVParameter();
483 const Standard_Real VMax = aBAS2.LastVParameter();
485 if( Precision::IsNegativeInfinite(VMin) ||
486 Precision::IsPositiveInfinite(VMax))
487 return Standard_True;
493 return Standard_True;
496 IntAna_QuadQuadGeo inter;
497 inter.Perform(aS1,aS2,Tolang,aTolSum, aHigh);
498 if(inter.TypeInter() == IntAna_Ellipse)
500 const gp_Elips anEl = inter.Ellipse(1);
501 const Standard_Real aMajorR = anEl.MajorRadius();
502 const Standard_Real aMinorR = anEl.MinorRadius();
504 return (aMajorR < 100000.0 * aMinorR);
508 return inter.IsDone();
511 //=======================================================================
513 //purpose : intersect surfaces of the faces
514 //=======================================================================
515 void IntTools_FaceFace::Perform(const TopoDS_Face& aF1,
516 const TopoDS_Face& aF2)
518 Standard_Boolean RestrictLine = Standard_False, hasCone = Standard_False;
520 if (myContext.IsNull()) {
521 myContext=new IntTools_Context;
524 mySeqOfCurve.Clear();
527 myIsDone = Standard_False;
533 const BRepAdaptor_Surface aBAS1(myFace1, Standard_False);
534 const BRepAdaptor_Surface aBAS2(myFace2, Standard_False);
535 GeomAbs_SurfaceType aType1=aBAS1.GetType();
536 GeomAbs_SurfaceType aType2=aBAS2.GetType();
538 const Standard_Boolean bReverse=SortTypes(aType1, aType2);
543 aType1=aBAS2.GetType();
544 aType2=aBAS1.GetType();
546 if (myListOfPnts.Extent())
548 Standard_Real aU1,aV1,aU2,aV2;
549 IntSurf_ListIteratorOfListOfPntOn2S aItP2S;
551 aItP2S.Initialize(myListOfPnts);
552 for (; aItP2S.More(); aItP2S.Next())
554 IntSurf_PntOn2S& aP2S=aItP2S.Value();
555 aP2S.Parameters(aU1,aV1,aU2,aV2);
556 aP2S.SetValue(aU2,aV2,aU1,aV1);
560 Standard_Boolean anAproxTmp = myApprox1;
561 myApprox1 = myApprox2;
562 myApprox2 = anAproxTmp;
566 const Handle(Geom_Surface) S1=BRep_Tool::Surface(myFace1);
567 const Handle(Geom_Surface) S2=BRep_Tool::Surface(myFace2);
569 const Standard_Real aTolF1=BRep_Tool::Tolerance(myFace1);
570 const Standard_Real aTolF2=BRep_Tool::Tolerance(myFace2);
572 Standard_Real TolArc = aTolF1 + aTolF2;
573 Standard_Real TolTang = TolArc;
575 const Standard_Boolean isFace1Quad = (aType1 == GeomAbs_Cylinder ||
576 aType1 == GeomAbs_Cone ||
577 aType1 == GeomAbs_Torus);
579 const Standard_Boolean isFace2Quad = (aType2 == GeomAbs_Cylinder ||
580 aType2 == GeomAbs_Cone ||
581 aType2 == GeomAbs_Torus);
583 if(aType1==GeomAbs_Plane && aType2==GeomAbs_Plane) {
584 Standard_Real umin, umax, vmin, vmax;
586 BRepTools::UVBounds(myFace1, umin, umax, vmin, vmax);
587 CorrectPlaneBoundaries(umin, umax, vmin, vmax);
588 myHS1->ChangeSurface().Load(S1, umin, umax, vmin, vmax);
590 BRepTools::UVBounds(myFace2, umin, umax, vmin, vmax);
591 CorrectPlaneBoundaries(umin, umax, vmin, vmax);
592 myHS2->ChangeSurface().Load(S2, umin, umax, vmin, vmax);
594 Standard_Real TolAng = 1.e-8;
596 PerformPlanes(myHS1, myHS2, TolAng, TolTang, myApprox1, myApprox2,
597 mySeqOfCurve, myTangentFaces);
599 myIsDone = Standard_True;
601 if(!myTangentFaces) {
602 const Standard_Integer NbLinPP = mySeqOfCurve.Length();
604 Standard_Real aTolFMax;
605 myTolReached3d = 1.e-7;
606 aTolFMax=Max(aTolF1, aTolF2);
607 if (aTolFMax>myTolReached3d) {
608 myTolReached3d=aTolFMax;
611 myTolReached2d = myTolReached3d;
614 Handle(Geom2d_Curve) aC2D1, aC2D2;
615 const Standard_Integer aNbLin = mySeqOfCurve.Length();
616 for (Standard_Integer i = 1; i <= aNbLin; ++i) {
617 IntTools_Curve& aIC=mySeqOfCurve(i);
618 aC2D1=aIC.FirstCurve2d();
619 aC2D2=aIC.SecondCurve2d();
620 aIC.SetFirstCurve2d(aC2D2);
621 aIC.SetSecondCurve2d(aC2D1);
627 }//if(aType1==GeomAbs_Plane && aType2==GeomAbs_Plane){
629 if ((aType1==GeomAbs_Plane) && isFace2Quad)
631 Standard_Real umin, umax, vmin, vmax;
633 BRepTools::UVBounds(myFace1, umin, umax, vmin, vmax);
634 CorrectPlaneBoundaries(umin, umax, vmin, vmax);
635 myHS1->ChangeSurface().Load(S1, umin, umax, vmin, vmax);
637 BRepTools::UVBounds(myFace2, umin, umax, vmin, vmax);
638 CorrectSurfaceBoundaries(myFace2, (aTolF1 + aTolF2) * 2., umin, umax, vmin, vmax);
639 myHS2->ChangeSurface().Load(S2, umin, umax, vmin, vmax);
641 if( aType2==GeomAbs_Cone ) {
643 hasCone = Standard_True;
646 else if ((aType2==GeomAbs_Plane) && isFace1Quad)
648 Standard_Real umin, umax, vmin, vmax;
650 BRepTools::UVBounds(myFace1, umin, umax, vmin, vmax);
651 CorrectSurfaceBoundaries(myFace1, (aTolF1 + aTolF2) * 2., umin, umax, vmin, vmax);
652 myHS1->ChangeSurface().Load(S1, umin, umax, vmin, vmax);
654 BRepTools::UVBounds(myFace2, umin, umax, vmin, vmax);
655 CorrectPlaneBoundaries(umin, umax, vmin, vmax);
656 myHS2->ChangeSurface().Load(S2, umin, umax, vmin, vmax);
658 if( aType1==GeomAbs_Cone ) {
660 hasCone = Standard_True;
665 Standard_Real umin, umax, vmin, vmax;
666 BRepTools::UVBounds(myFace1, umin, umax, vmin, vmax);
667 CorrectSurfaceBoundaries(myFace1, (aTolF1 + aTolF2) * 2., umin, umax, vmin, vmax);
668 myHS1->ChangeSurface().Load(S1, umin, umax, vmin, vmax);
669 BRepTools::UVBounds(myFace2, umin, umax, vmin, vmax);
670 CorrectSurfaceBoundaries(myFace2, (aTolF1 + aTolF2) * 2., umin, umax, vmin, vmax);
671 myHS2->ChangeSurface().Load(S2, umin, umax, vmin, vmax);
674 const Handle(IntTools_TopolTool) dom1 = new IntTools_TopolTool(myHS1);
675 const Handle(IntTools_TopolTool) dom2 = new IntTools_TopolTool(myHS2);
677 myLConstruct.Load(dom1, dom2, myHS1, myHS2);
680 Tolerances(myHS1, myHS2, TolTang);
683 const Standard_Real UVMaxStep = 0.001;
684 const Standard_Real Deflection = (hasCone) ? 0.085 : 0.1;
685 myIntersector.SetTolerances(TolArc, TolTang, UVMaxStep, Deflection);
688 if((myHS1->IsUClosed() && !myHS1->IsUPeriodic()) ||
689 (myHS1->IsVClosed() && !myHS1->IsVPeriodic()) ||
690 (myHS2->IsUClosed() && !myHS2->IsUPeriodic()) ||
691 (myHS2->IsVClosed() && !myHS2->IsVPeriodic()))
693 RestrictLine = Standard_True;
696 if((aType1 != GeomAbs_BSplineSurface) &&
697 (aType1 != GeomAbs_BezierSurface) &&
698 (aType1 != GeomAbs_OtherSurface) &&
699 (aType2 != GeomAbs_BSplineSurface) &&
700 (aType2 != GeomAbs_BezierSurface) &&
701 (aType2 != GeomAbs_OtherSurface))
703 RestrictLine = Standard_True;
705 if ((aType1 == GeomAbs_Torus) ||
706 (aType2 == GeomAbs_Torus))
708 myListOfPnts.Clear();
715 TopExp_Explorer aExp;
716 for(Standard_Integer i = 0; (!RestrictLine) && (i < 2); i++)
718 const TopoDS_Face& aF=(!i) ? myFace1 : myFace2;
719 aExp.Init(aF, TopAbs_EDGE);
720 for(; aExp.More(); aExp.Next())
722 const TopoDS_Edge& aE=TopoDS::Edge(aExp.Current());
724 if(BRep_Tool::Degenerated(aE))
726 RestrictLine = Standard_True;
733 const Standard_Boolean isGeomInt = isTreatAnalityc(aF1, aF2);
734 myIntersector.Perform(myHS1, dom1, myHS2, dom2, TolArc, TolTang,
735 myListOfPnts, RestrictLine, isGeomInt);
737 myIsDone = myIntersector.IsDone();
741 myTangentFaces=myIntersector.TangentFaces();
742 if (myTangentFaces) {
747 myListOfPnts.Clear(); // to use LineConstructor
750 const Standard_Integer aNbLinIntersector = myIntersector.NbLines();
751 for (Standard_Integer i=1; i <= aNbLinIntersector; ++i) {
752 MakeCurve(i, dom1, dom2);
755 ComputeTolReached3d();
758 Handle(Geom2d_Curve) aC2D1, aC2D2;
760 const Standard_Integer aNbLinSeqOfCurve =mySeqOfCurve.Length();
761 for (Standard_Integer i=1; i<=aNbLinSeqOfCurve; ++i)
763 IntTools_Curve& aIC=mySeqOfCurve(i);
764 aC2D1=aIC.FirstCurve2d();
765 aC2D2=aIC.SecondCurve2d();
766 aIC.SetFirstCurve2d(aC2D2);
767 aIC.SetSecondCurve2d(aC2D1);
772 Standard_Boolean bValid2D1, bValid2D2;
773 Standard_Real U1,V1,U2,V2;
774 IntTools_PntOnFace aPntOnF1, aPntOnF2;
775 IntTools_PntOn2Faces aPntOn2Faces;
777 const Standard_Integer aNbPnts = myIntersector.NbPnts();
778 for (Standard_Integer i=1; i <= aNbPnts; ++i)
780 const IntSurf_PntOn2S& aISPnt=myIntersector.Point(i).PntOn2S();
781 const gp_Pnt& aPnt=aISPnt.Value();
782 aISPnt.Parameters(U1,V1,U2,V2);
784 // check the validity of the intersection point for the faces
785 bValid2D1 = myContext->IsPointInOnFace(myFace1, gp_Pnt2d(U1, V1));
790 bValid2D2 = myContext->IsPointInOnFace(myFace2, gp_Pnt2d(U2, V2));
795 // add the intersection point
796 aPntOnF1.Init(myFace1, aPnt, U1, V1);
797 aPntOnF2.Init(myFace2, aPnt, U2, V2);
801 aPntOn2Faces.SetP1(aPntOnF1);
802 aPntOn2Faces.SetP2(aPntOnF2);
806 aPntOn2Faces.SetP2(aPntOnF1);
807 aPntOn2Faces.SetP1(aPntOnF2);
810 myPnts.Append(aPntOn2Faces);
815 //=======================================================================
816 //function : ComputeTolerance
818 //=======================================================================
819 Standard_Real IntTools_FaceFace::ComputeTolerance()
821 Standard_Integer i, j, aNbLin;
822 Standard_Real aFirst, aLast, aD, aDMax, aT;
823 Handle(Geom_Surface) aS1, aS2;
826 aNbLin = mySeqOfCurve.Length();
828 aS1 = myHS1->ChangeSurface().Surface();
829 aS2 = myHS2->ChangeSurface().Surface();
831 for (i = 1; i <= aNbLin; ++i)
833 const IntTools_Curve& aIC = mySeqOfCurve(i);
834 const Handle(Geom_Curve)& aC3D = aIC.Curve();
840 aFirst = aC3D->FirstParameter();
841 aLast = aC3D->LastParameter();
843 const Handle(Geom2d_Curve)& aC2D1 = aIC.FirstCurve2d();
844 const Handle(Geom2d_Curve)& aC2D2 = aIC.SecondCurve2d();
846 for (j = 0; j < 2; ++j)
848 const Handle(Geom2d_Curve)& aC2D = !j ? aC2D1 : aC2D2;
849 const Handle(Geom_Surface)& aS = !j ? aS1 : aS2;
853 if (IntTools_Tools::ComputeTolerance
854 (aC3D, aC2D, aS, aFirst, aLast, aD, aT))
864 const TopoDS_Face& aF = !j ? myFace1 : myFace2;
865 aD = FindMaxDistance(aC3D, aFirst, aLast, aF, myContext);
877 //=======================================================================
878 //function :ComputeTolReached3d
880 //=======================================================================
881 void IntTools_FaceFace::ComputeTolReached3d()
883 Standard_Integer aNbLin;
884 GeomAbs_SurfaceType aType1, aType2;
886 aNbLin=myIntersector.NbLines();
891 aType1=myHS1->Surface().GetType();
892 aType2=myHS2->Surface().GetType();
894 if (aType1==GeomAbs_Cylinder && aType2==GeomAbs_Cylinder)
898 Handle(IntPatch_Line) aIL1, aIL2;
899 IntPatch_IType aTL1, aTL2;
901 aIL1=myIntersector.Line(1);
902 aIL2=myIntersector.Line(2);
903 aTL1=aIL1->ArcType();
904 aTL2=aIL2->ArcType();
905 if (aTL1==IntPatch_Lin && aTL2==IntPatch_Lin) {
906 Standard_Real aD, aDTresh, dTol;
912 aL1=Handle(IntPatch_GLine)::DownCast(aIL1)->Line();
913 aL2=Handle(IntPatch_GLine)::DownCast(aIL2)->Line();
914 aD=aL1.Distance(aL2);
917 {//In order to avoid creation too thin face
918 myTolReached3d=aD+dTol;
922 }// if (aType1==GeomAbs_Cylinder && aType2==GeomAbs_Cylinder) {
925 Standard_Real aDMax = ComputeTolerance();
926 if (aDMax > myTolReached3d)
928 myTolReached3d = aDMax;
932 //=======================================================================
933 //function : MakeCurve
935 //=======================================================================
936 void IntTools_FaceFace::MakeCurve(const Standard_Integer Index,
937 const Handle(Adaptor3d_TopolTool)& dom1,
938 const Handle(Adaptor3d_TopolTool)& dom2)
940 Standard_Boolean bDone, rejectSurface, reApprox, bAvoidLineConstructor;
941 Standard_Boolean ok, bPCurvesOk;
942 Standard_Integer i, j, aNbParts;
943 Standard_Real fprm, lprm;
945 Handle(IntPatch_Line) L;
947 Handle(Geom_Curve) newc;
949 const Standard_Real TOLCHECK =0.0000001;
950 const Standard_Real TOLANGCHECK=0.1;
952 rejectSurface = Standard_False;
953 reApprox = Standard_False;
955 bPCurvesOk = Standard_True;
960 bAvoidLineConstructor = Standard_False;
961 L = myIntersector.Line(Index);
964 if(typl==IntPatch_Walking) {
965 Handle(IntPatch_Line) anewL;
967 Handle(IntPatch_WLine) aWLine (Handle(IntPatch_WLine)::DownCast(L));
968 anewL = ComputePurgedWLine(aWLine, myHS1, myHS2);
974 //Handle(IntPatch_WLine) aWLineX (Handle(IntPatch_WLine)::DownCast(L));
975 //DumpWLine(aWLineX);
978 if(!myListOfPnts.IsEmpty()) {
979 bAvoidLineConstructor = Standard_True;
982 Standard_Integer nbp = aWLine->NbPnts();
983 const IntSurf_PntOn2S& p1 = aWLine->Point(1);
984 const IntSurf_PntOn2S& p2 = aWLine->Point(nbp);
986 const gp_Pnt& P1 = p1.Value();
987 const gp_Pnt& P2 = p2.Value();
989 if(P1.SquareDistance(P2) < 1.e-14) {
990 bAvoidLineConstructor = Standard_False;
998 if(!bAvoidLineConstructor) {
999 myLConstruct.Perform(L);
1001 bDone=myLConstruct.IsDone();
1007 if(typl != IntPatch_Restriction)
1009 aNbParts=myLConstruct.NbParts();
1020 //########################################
1021 // Line, Parabola, Hyperbola
1022 //########################################
1024 case IntPatch_Parabola:
1025 case IntPatch_Hyperbola: {
1026 if (typl == IntPatch_Lin) {
1028 new Geom_Line (Handle(IntPatch_GLine)::DownCast(L)->Line());
1031 else if (typl == IntPatch_Parabola) {
1033 new Geom_Parabola(Handle(IntPatch_GLine)::DownCast(L)->Parabola());
1036 else if (typl == IntPatch_Hyperbola) {
1038 new Geom_Hyperbola (Handle(IntPatch_GLine)::DownCast(L)->Hyperbola());
1042 if (typl == IntPatch_Lin) {
1043 TolR3d (myFace1, myFace2, myTolReached3d);
1046 aNbParts=myLConstruct.NbParts();
1047 for (i=1; i<=aNbParts; i++) {
1048 Standard_Boolean bFNIt, bLPIt;
1050 myLConstruct.Part(i, fprm, lprm);
1052 bFNIt=Precision::IsNegativeInfinite(fprm);
1053 bLPIt=Precision::IsPositiveInfinite(lprm);
1055 if (!bFNIt && !bLPIt) {
1057 IntTools_Curve aCurve;
1059 Handle(Geom_TrimmedCurve) aCT3D=new Geom_TrimmedCurve(newc, fprm, lprm);
1060 aCurve.SetCurve(aCT3D);
1061 if (typl == IntPatch_Parabola) {
1062 Standard_Real aTolF1, aTolF2, aTolBase;
1064 aTolF1 = BRep_Tool::Tolerance(myFace1);
1065 aTolF2 = BRep_Tool::Tolerance(myFace2);
1066 aTolBase=aTolF1+aTolF2;
1067 myTolReached3d=IntTools_Tools::CurveTolerance(aCT3D, aTolBase);
1070 aCurve.SetCurve(new Geom_TrimmedCurve(newc, fprm, lprm));
1072 Handle (Geom2d_Curve) C2d;
1073 GeomInt_IntSS::BuildPCurves(fprm, lprm, Tolpc,
1074 myHS1->ChangeSurface().Surface(), newc, C2d);
1075 if(Tolpc>myTolReached2d || myTolReached2d==0.) {
1076 myTolReached2d=Tolpc;
1079 aCurve.SetFirstCurve2d(new Geom2d_TrimmedCurve(C2d,fprm,lprm));
1082 Handle(Geom2d_BSplineCurve) H1;
1084 aCurve.SetFirstCurve2d(H1);
1088 Handle (Geom2d_Curve) C2d;
1089 GeomInt_IntSS::BuildPCurves(fprm, lprm, Tolpc,
1090 myHS2->ChangeSurface().Surface(), newc, C2d);
1091 if(Tolpc>myTolReached2d || myTolReached2d==0.) {
1092 myTolReached2d=Tolpc;
1095 aCurve.SetSecondCurve2d(new Geom2d_TrimmedCurve(C2d,fprm,lprm));
1098 Handle(Geom2d_BSplineCurve) H1;
1100 aCurve.SetSecondCurve2d(H1);
1102 mySeqOfCurve.Append(aCurve);
1103 } //if (!bFNIt && !bLPIt) {
1105 // on regarde si on garde
1107 Standard_Real aTestPrm, dT=100.;
1110 if (bFNIt && !bLPIt) {
1113 else if (!bFNIt && bLPIt) {
1117 // i.e, if (bFNIt && bLPIt)
1118 aTestPrm=IntTools_Tools::IntermediatePoint(-dT, dT);
1121 gp_Pnt ptref(newc->Value(aTestPrm));
1123 GeomAbs_SurfaceType typS1 = myHS1->GetType();
1124 GeomAbs_SurfaceType typS2 = myHS2->GetType();
1125 if( typS1 == GeomAbs_SurfaceOfExtrusion ||
1126 typS1 == GeomAbs_OffsetSurface ||
1127 typS1 == GeomAbs_SurfaceOfRevolution ||
1128 typS2 == GeomAbs_SurfaceOfExtrusion ||
1129 typS2 == GeomAbs_OffsetSurface ||
1130 typS2 == GeomAbs_SurfaceOfRevolution) {
1131 Handle(Geom2d_BSplineCurve) H1;
1132 mySeqOfCurve.Append(IntTools_Curve(newc, H1, H1));
1136 Standard_Real u1, v1, u2, v2, Tol;
1138 Tol = Precision::Confusion();
1139 Parameters(myHS1, myHS2, ptref, u1, v1, u2, v2);
1140 ok = (dom1->Classify(gp_Pnt2d(u1, v1), Tol) != TopAbs_OUT);
1142 ok = (dom2->Classify(gp_Pnt2d(u2,v2),Tol) != TopAbs_OUT);
1145 Handle(Geom2d_BSplineCurve) H1;
1146 mySeqOfCurve.Append(IntTools_Curve(newc, H1, H1));
1149 }// for (i=1; i<=aNbParts; i++) {
1150 }// case IntPatch_Lin: case IntPatch_Parabola: case IntPatch_Hyperbola:
1153 //########################################
1154 // Circle and Ellipse
1155 //########################################
1156 case IntPatch_Circle:
1157 case IntPatch_Ellipse: {
1159 if (typl == IntPatch_Circle) {
1160 newc = new Geom_Circle
1161 (Handle(IntPatch_GLine)::DownCast(L)->Circle());
1163 else { //IntPatch_Ellipse
1164 newc = new Geom_Ellipse
1165 (Handle(IntPatch_GLine)::DownCast(L)->Ellipse());
1169 TolR3d (myFace1, myFace2, myTolReached3d);
1171 aNbParts=myLConstruct.NbParts();
1173 Standard_Real aPeriod, aNul;
1174 TColStd_SequenceOfReal aSeqFprm, aSeqLprm;
1179 for (i=1; i<=aNbParts; i++) {
1180 myLConstruct.Part(i, fprm, lprm);
1182 if (fprm < aNul && lprm > aNul) {
1183 // interval that goes through 0. is divided on two intervals;
1184 while (fprm<aNul || fprm>aPeriod) fprm=fprm+aPeriod;
1185 while (lprm<aNul || lprm>aPeriod) lprm=lprm+aPeriod;
1187 if((aPeriod - fprm) > Tolpc) {
1188 aSeqFprm.Append(fprm);
1189 aSeqLprm.Append(aPeriod);
1192 gp_Pnt P1 = newc->Value(fprm);
1193 gp_Pnt P2 = newc->Value(aPeriod);
1194 Standard_Real aTolDist = BRep_Tool::Tolerance(myFace1) + BRep_Tool::Tolerance(myFace2);
1195 aTolDist = (myTolReached3d > aTolDist) ? myTolReached3d : aTolDist;
1197 if(P1.Distance(P2) > aTolDist) {
1198 Standard_Real anewpar = fprm;
1200 if(ParameterOutOfBoundary(fprm, newc, myFace1, myFace2,
1201 lprm, Standard_False, anewpar, myContext)) {
1204 aSeqFprm.Append(fprm);
1205 aSeqLprm.Append(aPeriod);
1210 if((lprm - aNul) > Tolpc) {
1211 aSeqFprm.Append(aNul);
1212 aSeqLprm.Append(lprm);
1215 gp_Pnt P1 = newc->Value(aNul);
1216 gp_Pnt P2 = newc->Value(lprm);
1217 Standard_Real aTolDist = BRep_Tool::Tolerance(myFace1) + BRep_Tool::Tolerance(myFace2);
1218 aTolDist = (myTolReached3d > aTolDist) ? myTolReached3d : aTolDist;
1220 if(P1.Distance(P2) > aTolDist) {
1221 Standard_Real anewpar = lprm;
1223 if(ParameterOutOfBoundary(lprm, newc, myFace1, myFace2,
1224 fprm, Standard_True, anewpar, myContext)) {
1227 aSeqFprm.Append(aNul);
1228 aSeqLprm.Append(lprm);
1234 aSeqFprm.Append(fprm);
1235 aSeqLprm.Append(lprm);
1240 aNbParts=aSeqFprm.Length();
1241 for (i=1; i<=aNbParts; i++) {
1245 Standard_Real aRealEpsilon=RealEpsilon();
1246 if (Abs(fprm) > aRealEpsilon || Abs(lprm-2.*M_PI) > aRealEpsilon) {
1247 //==============================================
1249 IntTools_Curve aCurve;
1250 Handle(Geom_TrimmedCurve) aTC3D=new Geom_TrimmedCurve(newc,fprm,lprm);
1251 aCurve.SetCurve(aTC3D);
1252 fprm=aTC3D->FirstParameter();
1253 lprm=aTC3D->LastParameter ();
1255 if (typl == IntPatch_Circle || typl == IntPatch_Ellipse) {////
1257 Handle (Geom2d_Curve) C2d;
1258 GeomInt_IntSS::BuildPCurves(fprm, lprm, Tolpc,
1259 myHS1->ChangeSurface().Surface(), newc, C2d);
1260 if(Tolpc>myTolReached2d || myTolReached2d==0) {
1261 myTolReached2d=Tolpc;
1264 aCurve.SetFirstCurve2d(C2d);
1267 Handle(Geom2d_BSplineCurve) H1;
1268 aCurve.SetFirstCurve2d(H1);
1273 Handle (Geom2d_Curve) C2d;
1274 GeomInt_IntSS::BuildPCurves(fprm,lprm,Tolpc,
1275 myHS2->ChangeSurface().Surface(),newc,C2d);
1276 if(Tolpc>myTolReached2d || myTolReached2d==0) {
1277 myTolReached2d=Tolpc;
1280 aCurve.SetSecondCurve2d(C2d);
1283 Handle(Geom2d_BSplineCurve) H1;
1284 aCurve.SetSecondCurve2d(H1);
1289 Handle(Geom2d_BSplineCurve) H1;
1290 aCurve.SetFirstCurve2d(H1);
1291 aCurve.SetSecondCurve2d(H1);
1293 mySeqOfCurve.Append(aCurve);
1294 //==============================================
1295 } //if (Abs(fprm) > RealEpsilon() || Abs(lprm-2.*M_PI) > RealEpsilon())
1298 // on regarde si on garde
1301 // if (Abs(fprm) < RealEpsilon() && Abs(lprm-2.*M_PI) < RealEpsilon()) {
1302 if (Abs(fprm) <= aRealEpsilon && Abs(lprm-2.*M_PI) <= aRealEpsilon) {
1303 IntTools_Curve aCurve;
1304 Handle(Geom_TrimmedCurve) aTC3D=new Geom_TrimmedCurve(newc,fprm,lprm);
1305 aCurve.SetCurve(aTC3D);
1306 fprm=aTC3D->FirstParameter();
1307 lprm=aTC3D->LastParameter ();
1310 Handle (Geom2d_Curve) C2d;
1311 GeomInt_IntSS::BuildPCurves(fprm,lprm,Tolpc,
1312 myHS1->ChangeSurface().Surface(),newc,C2d);
1313 if(Tolpc>myTolReached2d || myTolReached2d==0) {
1314 myTolReached2d=Tolpc;
1317 aCurve.SetFirstCurve2d(C2d);
1320 Handle(Geom2d_BSplineCurve) H1;
1321 aCurve.SetFirstCurve2d(H1);
1325 Handle (Geom2d_Curve) C2d;
1326 GeomInt_IntSS::BuildPCurves(fprm,lprm,Tolpc,
1327 myHS2->ChangeSurface().Surface(),newc,C2d);
1328 if(Tolpc>myTolReached2d || myTolReached2d==0) {
1329 myTolReached2d=Tolpc;
1332 aCurve.SetSecondCurve2d(C2d);
1335 Handle(Geom2d_BSplineCurve) H1;
1336 aCurve.SetSecondCurve2d(H1);
1338 mySeqOfCurve.Append(aCurve);
1343 Standard_Real aTwoPIdiv17, u1, v1, u2, v2, Tol;
1345 aTwoPIdiv17=2.*M_PI/17.;
1347 for (j=0; j<=17; j++) {
1348 gp_Pnt ptref (newc->Value (j*aTwoPIdiv17));
1349 Tol = Precision::Confusion();
1351 Parameters(myHS1, myHS2, ptref, u1, v1, u2, v2);
1352 ok = (dom1->Classify(gp_Pnt2d(u1,v1),Tol) != TopAbs_OUT);
1354 ok = (dom2->Classify(gp_Pnt2d(u2,v2),Tol) != TopAbs_OUT);
1357 IntTools_Curve aCurve;
1358 aCurve.SetCurve(newc);
1359 //==============================================
1360 if (typl == IntPatch_Circle || typl == IntPatch_Ellipse) {
1363 Handle (Geom2d_Curve) C2d;
1364 GeomInt_IntSS::BuildPCurves(fprm, lprm, Tolpc,
1365 myHS1->ChangeSurface().Surface(), newc, C2d);
1366 if(Tolpc>myTolReached2d || myTolReached2d==0) {
1367 myTolReached2d=Tolpc;
1370 aCurve.SetFirstCurve2d(C2d);
1373 Handle(Geom2d_BSplineCurve) H1;
1374 aCurve.SetFirstCurve2d(H1);
1378 Handle (Geom2d_Curve) C2d;
1379 GeomInt_IntSS::BuildPCurves(fprm, lprm, Tolpc,
1380 myHS2->ChangeSurface().Surface(), newc, C2d);
1381 if(Tolpc>myTolReached2d || myTolReached2d==0) {
1382 myTolReached2d=Tolpc;
1385 aCurve.SetSecondCurve2d(C2d);
1389 Handle(Geom2d_BSplineCurve) H1;
1390 aCurve.SetSecondCurve2d(H1);
1392 }// end of if (typl == IntPatch_Circle || typl == IntPatch_Ellipse)
1395 Handle(Geom2d_BSplineCurve) H1;
1397 aCurve.SetFirstCurve2d(H1);
1398 aCurve.SetSecondCurve2d(H1);
1400 //==============================================
1402 mySeqOfCurve.Append(aCurve);
1405 }// end of if (ok) {
1406 }// end of for (Standard_Integer j=0; j<=17; j++)
1407 }// end of else { on regarde si on garde
1408 }// for (i=1; i<=myLConstruct.NbParts(); i++)
1409 }// IntPatch_Circle: IntPatch_Ellipse:
1412 case IntPatch_Analytic: {
1413 IntSurf_Quadric quad1,quad2;
1414 GeomAbs_SurfaceType typs = myHS1->Surface().GetType();
1418 quad1.SetValue(myHS1->Surface().Plane());
1420 case GeomAbs_Cylinder:
1421 quad1.SetValue(myHS1->Surface().Cylinder());
1424 quad1.SetValue(myHS1->Surface().Cone());
1426 case GeomAbs_Sphere:
1427 quad1.SetValue(myHS1->Surface().Sphere());
1430 quad1.SetValue(myHS1->Surface().Torus());
1433 Standard_ConstructionError::Raise("GeomInt_IntSS::MakeCurve 1");
1436 typs = myHS2->Surface().GetType();
1440 quad2.SetValue(myHS2->Surface().Plane());
1442 case GeomAbs_Cylinder:
1443 quad2.SetValue(myHS2->Surface().Cylinder());
1446 quad2.SetValue(myHS2->Surface().Cone());
1448 case GeomAbs_Sphere:
1449 quad2.SetValue(myHS2->Surface().Sphere());
1452 quad2.SetValue(myHS2->Surface().Torus());
1455 Standard_ConstructionError::Raise("GeomInt_IntSS::MakeCurve 2");
1459 IntPatch_ALineToWLine convert (quad1, quad2);
1462 aNbParts=myLConstruct.NbParts();
1463 for (i=1; i<=aNbParts; i++) {
1464 myLConstruct.Part(i, fprm, lprm);
1465 Handle(IntPatch_WLine) WL =
1466 convert.MakeWLine(Handle(IntPatch_ALine)::DownCast(L), fprm, lprm);
1468 Handle(Geom2d_BSplineCurve) H1;
1469 Handle(Geom2d_BSplineCurve) H2;
1472 H1 = MakeBSpline2d(WL, 1, WL->NbPnts(), Standard_True);
1476 H2 = MakeBSpline2d(WL, 1, WL->NbPnts(), Standard_False);
1479 mySeqOfCurve.Append(IntTools_Curve(MakeBSpline(WL,1,WL->NbPnts()), H1, H2));
1483 else { // myApprox=TRUE
1484 GeomInt_WLApprox theapp3d;
1486 Standard_Real tol2d = myTolApprox;
1488 theapp3d.SetParameters(myTolApprox, tol2d, 4, 8, 0, Standard_True);
1490 aNbParts=myLConstruct.NbParts();
1491 for (i=1; i<=aNbParts; i++) {
1492 myLConstruct.Part(i, fprm, lprm);
1493 Handle(IntPatch_WLine) WL =
1494 convert.MakeWLine(Handle(IntPatch_ALine):: DownCast(L),fprm,lprm);
1496 theapp3d.Perform(myHS1,myHS2,WL,Standard_True,myApprox1,myApprox2, 1, WL->NbPnts());
1498 if (!theapp3d.IsDone()) {
1500 Handle(Geom2d_BSplineCurve) H1;
1501 Handle(Geom2d_BSplineCurve) H2;
1504 H1 = MakeBSpline2d(WL, 1, WL->NbPnts(), Standard_True);
1508 H2 = MakeBSpline2d(WL, 1, WL->NbPnts(), Standard_False);
1511 mySeqOfCurve.Append(IntTools_Curve(MakeBSpline(WL,1,WL->NbPnts()), H1, H2));
1515 if(myApprox1 || myApprox2) {
1516 if( theapp3d.TolReached2d()>myTolReached2d || myTolReached2d==0) {
1517 myTolReached2d = theapp3d.TolReached2d();
1521 if( theapp3d.TolReached3d()>myTolReached3d || myTolReached3d==0) {
1522 myTolReached3d = theapp3d.TolReached3d();
1525 Standard_Integer aNbMultiCurves, nbpoles;
1526 aNbMultiCurves=theapp3d.NbMultiCurves();
1527 for (j=1; j<=aNbMultiCurves; j++) {
1528 const AppParCurves_MultiBSpCurve& mbspc = theapp3d.Value(j);
1529 nbpoles = mbspc.NbPoles();
1531 TColgp_Array1OfPnt tpoles(1, nbpoles);
1532 mbspc.Curve(1, tpoles);
1533 Handle(Geom_BSplineCurve) BS=new Geom_BSplineCurve(tpoles,
1535 mbspc.Multiplicities(),
1538 GeomLib_CheckBSplineCurve Check(BS,TOLCHECK,TOLANGCHECK);
1539 Check.FixTangent(Standard_True,Standard_True);
1541 IntTools_Curve aCurve;
1542 aCurve.SetCurve(BS);
1545 TColgp_Array1OfPnt2d tpoles2d(1,nbpoles);
1546 mbspc.Curve(2,tpoles2d);
1547 Handle(Geom2d_BSplineCurve) BS2=new Geom2d_BSplineCurve(tpoles2d,
1549 mbspc.Multiplicities(),
1552 GeomLib_Check2dBSplineCurve newCheck(BS2,TOLCHECK,TOLANGCHECK);
1553 newCheck.FixTangent(Standard_True,Standard_True);
1555 aCurve.SetFirstCurve2d(BS2);
1558 Handle(Geom2d_BSplineCurve) H1;
1559 aCurve.SetFirstCurve2d(H1);
1563 TColgp_Array1OfPnt2d tpoles2d(1, nbpoles);
1564 Standard_Integer TwoOrThree;
1565 TwoOrThree=myApprox1 ? 3 : 2;
1566 mbspc.Curve(TwoOrThree, tpoles2d);
1567 Handle(Geom2d_BSplineCurve) BS2 =new Geom2d_BSplineCurve(tpoles2d,
1569 mbspc.Multiplicities(),
1572 GeomLib_Check2dBSplineCurve newCheck(BS2,TOLCHECK,TOLANGCHECK);
1573 newCheck.FixTangent(Standard_True,Standard_True);
1575 aCurve.SetSecondCurve2d(BS2);
1578 Handle(Geom2d_BSplineCurve) H2;
1579 aCurve.SetSecondCurve2d(H2);
1582 mySeqOfCurve.Append(aCurve);
1584 }// for (j=1; j<=aNbMultiCurves; j++) {
1585 }// else from if (!theapp3d.IsDone())
1586 }// for (i=1; i<=aNbParts; i++) {
1587 }// else { // myApprox=TRUE
1588 }// case IntPatch_Analytic:
1591 case IntPatch_Walking:{
1592 Handle(IntPatch_WLine) WL =
1593 Handle(IntPatch_WLine)::DownCast(L);
1600 Standard_Integer ifprm, ilprm;
1604 if(!bAvoidLineConstructor){
1605 aNbParts=myLConstruct.NbParts();
1607 for (i=1; i<=aNbParts; ++i) {
1608 Handle(Geom2d_BSplineCurve) H1, H2;
1609 Handle(Geom_Curve) aBSp;
1611 if(bAvoidLineConstructor) {
1613 ilprm = WL->NbPnts();
1616 myLConstruct.Part(i, fprm, lprm);
1617 ifprm=(Standard_Integer)fprm;
1618 ilprm=(Standard_Integer)lprm;
1622 H1 = MakeBSpline2d(WL, ifprm, ilprm, Standard_True);
1626 H2 = MakeBSpline2d(WL, ifprm, ilprm, Standard_False);
1629 aBSp=MakeBSpline(WL, ifprm, ilprm);
1630 IntTools_Curve aIC(aBSp, H1, H2);
1631 mySeqOfCurve.Append(aIC);
1632 }// for (i=1; i<=aNbParts; ++i) {
1633 }// if (!myApprox) {
1636 Standard_Boolean bIsDecomposited;
1637 Standard_Integer nbiter, aNbSeqOfL;
1638 Standard_Real tol2d, aTolApproxImp;
1639 IntPatch_SequenceOfLine aSeqOfL;
1640 GeomInt_WLApprox theapp3d;
1641 Approx_ParametrizationType aParType = Approx_ChordLength;
1643 Standard_Boolean anApprox1 = myApprox1;
1644 Standard_Boolean anApprox2 = myApprox2;
1646 aTolApproxImp=1.e-5;
1647 tol2d = myTolApprox;
1649 GeomAbs_SurfaceType typs1, typs2;
1650 typs1 = myHS1->Surface().GetType();
1651 typs2 = myHS2->Surface().GetType();
1652 Standard_Boolean anWithPC = Standard_True;
1654 if(typs1 == GeomAbs_Cylinder && typs2 == GeomAbs_Sphere) {
1656 ApproxWithPCurves(myHS1->Surface().Cylinder(), myHS2->Surface().Sphere());
1658 else if (typs1 == GeomAbs_Sphere && typs2 == GeomAbs_Cylinder) {
1660 ApproxWithPCurves(myHS2->Surface().Cylinder(), myHS1->Surface().Sphere());
1664 myTolApprox = aTolApproxImp;//1.e-5;
1665 anApprox1 = Standard_False;
1666 anApprox2 = Standard_False;
1668 tol2d = myTolApprox;
1671 if(myHS1 == myHS2) {
1672 theapp3d.SetParameters(myTolApprox, tol2d, 4, 8, 0, Standard_False, aParType);
1673 rejectSurface = Standard_True;
1676 if(reApprox && !rejectSurface)
1677 theapp3d.SetParameters(myTolApprox, tol2d, 4, 8, 0, Standard_False, aParType);
1679 Standard_Integer iDegMax, iDegMin, iNbIter;
1681 ApproxParameters(myHS1, myHS2, iDegMin, iDegMax, iNbIter);
1682 theapp3d.SetParameters(myTolApprox, tol2d, iDegMin, iDegMax, iNbIter, Standard_True, aParType);
1686 Standard_Real aReachedTol = Precision::Confusion();
1687 bIsDecomposited=DecompositionOfWLine(WL,
1693 bAvoidLineConstructor,
1697 if ( bIsDecomposited && ( myTolReached3d < aReachedTol ) ) {
1698 myTolReached3d = aReachedTol;
1701 aNbSeqOfL=aSeqOfL.Length();
1703 if (bIsDecomposited) {
1709 if (!bAvoidLineConstructor) {
1710 aNbParts=myLConstruct.NbParts();
1715 for(i = 1; i <= nbiter; ++i) {
1716 if(bIsDecomposited) {
1717 WL = Handle(IntPatch_WLine)::DownCast(aSeqOfL.Value(i));
1719 ilprm = WL->NbPnts();
1722 if(bAvoidLineConstructor) {
1724 ilprm = WL->NbPnts();
1727 myLConstruct.Part(i, fprm, lprm);
1728 ifprm = (Standard_Integer)fprm;
1729 ilprm = (Standard_Integer)lprm;
1733 //-- Si une des surfaces est un plan , on approxime en 2d
1734 //-- sur cette surface et on remonte les points 2d en 3d.
1735 if(typs1 == GeomAbs_Plane) {
1736 theapp3d.Perform(myHS1, myHS2, WL, Standard_False,Standard_True, myApprox2,ifprm,ilprm);
1738 else if(typs2 == GeomAbs_Plane) {
1739 theapp3d.Perform(myHS1,myHS2,WL,Standard_False,myApprox1,Standard_True,ifprm,ilprm);
1743 if (myHS1 != myHS2){
1744 if ((typs1==GeomAbs_BezierSurface || typs1==GeomAbs_BSplineSurface) &&
1745 (typs2==GeomAbs_BezierSurface || typs2==GeomAbs_BSplineSurface)) {
1747 theapp3d.SetParameters(myTolApprox, tol2d, 4, 8, 0, Standard_True, aParType);
1749 Standard_Boolean bUseSurfaces;
1750 bUseSurfaces=NotUseSurfacesForApprox(myFace1, myFace2, WL, ifprm, ilprm);
1753 rejectSurface = Standard_True;
1755 theapp3d.SetParameters(myTolApprox, tol2d, 4, 8, 0, Standard_False, aParType);
1760 theapp3d.Perform(myHS1,myHS2,WL,Standard_True,anApprox1,anApprox2,ifprm,ilprm);
1763 if (!theapp3d.IsDone()) {
1764 Handle(Geom2d_BSplineCurve) H1;
1765 Handle(Geom2d_BSplineCurve) H2;
1767 Handle(Geom_Curve) aBSp=MakeBSpline(WL,ifprm, ilprm);
1770 H1 = MakeBSpline2d(WL, ifprm, ilprm, Standard_True);
1774 H2 = MakeBSpline2d(WL, ifprm, ilprm, Standard_False);
1777 IntTools_Curve aIC(aBSp, H1, H2);
1778 mySeqOfCurve.Append(aIC);
1782 if(myApprox1 || myApprox2 || (typs1==GeomAbs_Plane || typs2==GeomAbs_Plane)) {
1783 if( theapp3d.TolReached2d()>myTolReached2d || myTolReached2d==0.) {
1784 myTolReached2d = theapp3d.TolReached2d();
1787 if(typs1==GeomAbs_Plane || typs2==GeomAbs_Plane) {
1788 myTolReached3d = myTolReached2d;
1790 if (typs1==GeomAbs_Torus || typs2==GeomAbs_Torus) {
1791 if (myTolReached3d<1.e-6) {
1792 myTolReached3d = theapp3d.TolReached3d();
1793 myTolReached3d=1.e-6;
1797 else if( theapp3d.TolReached3d()>myTolReached3d || myTolReached3d==0.) {
1798 myTolReached3d = theapp3d.TolReached3d();
1801 Standard_Integer aNbMultiCurves, nbpoles;
1802 aNbMultiCurves=theapp3d.NbMultiCurves();
1803 for (j=1; j<=aNbMultiCurves; j++) {
1804 if(typs1 == GeomAbs_Plane) {
1805 const AppParCurves_MultiBSpCurve& mbspc = theapp3d.Value(j);
1806 nbpoles = mbspc.NbPoles();
1808 TColgp_Array1OfPnt2d tpoles2d(1,nbpoles);
1809 TColgp_Array1OfPnt tpoles(1,nbpoles);
1811 mbspc.Curve(1,tpoles2d);
1812 const gp_Pln& Pln = myHS1->Surface().Plane();
1814 Standard_Integer ik;
1815 for(ik = 1; ik<= nbpoles; ik++) {
1817 ElSLib::Value(tpoles2d.Value(ik).X(),
1818 tpoles2d.Value(ik).Y(),
1822 Handle(Geom_BSplineCurve) BS =
1823 new Geom_BSplineCurve(tpoles,
1825 mbspc.Multiplicities(),
1827 GeomLib_CheckBSplineCurve Check(BS,TOLCHECK,TOLANGCHECK);
1828 Check.FixTangent(Standard_True, Standard_True);
1830 IntTools_Curve aCurve;
1831 aCurve.SetCurve(BS);
1834 Handle(Geom2d_BSplineCurve) BS1 =
1835 new Geom2d_BSplineCurve(tpoles2d,
1837 mbspc.Multiplicities(),
1839 GeomLib_Check2dBSplineCurve Check1(BS1,TOLCHECK,TOLANGCHECK);
1840 Check1.FixTangent(Standard_True,Standard_True);
1842 // ############################################
1843 if(!rejectSurface && !reApprox) {
1844 Standard_Boolean isValid = IsCurveValid(BS1);
1846 reApprox = Standard_True;
1850 // ############################################
1851 aCurve.SetFirstCurve2d(BS1);
1854 Handle(Geom2d_BSplineCurve) H1;
1855 aCurve.SetFirstCurve2d(H1);
1859 mbspc.Curve(2, tpoles2d);
1861 Handle(Geom2d_BSplineCurve) BS2 = new Geom2d_BSplineCurve(tpoles2d,
1863 mbspc.Multiplicities(),
1865 GeomLib_Check2dBSplineCurve newCheck(BS2,TOLCHECK,TOLANGCHECK);
1866 newCheck.FixTangent(Standard_True,Standard_True);
1868 // ###########################################
1869 if(!rejectSurface && !reApprox) {
1870 Standard_Boolean isValid = IsCurveValid(BS2);
1872 reApprox = Standard_True;
1876 // ###########################################
1878 aCurve.SetSecondCurve2d(BS2);
1881 Handle(Geom2d_BSplineCurve) H2;
1883 aCurve.SetSecondCurve2d(H2);
1886 mySeqOfCurve.Append(aCurve);
1888 }//if(typs1 == GeomAbs_Plane) {
1890 else if(typs2 == GeomAbs_Plane)
1892 const AppParCurves_MultiBSpCurve& mbspc = theapp3d.Value(j);
1893 nbpoles = mbspc.NbPoles();
1895 TColgp_Array1OfPnt2d tpoles2d(1,nbpoles);
1896 TColgp_Array1OfPnt tpoles(1,nbpoles);
1897 mbspc.Curve((myApprox1==Standard_True)? 2 : 1,tpoles2d);
1898 const gp_Pln& Pln = myHS2->Surface().Plane();
1900 Standard_Integer ik;
1901 for(ik = 1; ik<= nbpoles; ik++) {
1903 ElSLib::Value(tpoles2d.Value(ik).X(),
1904 tpoles2d.Value(ik).Y(),
1909 Handle(Geom_BSplineCurve) BS=new Geom_BSplineCurve(tpoles,
1911 mbspc.Multiplicities(),
1913 GeomLib_CheckBSplineCurve Check(BS,TOLCHECK,TOLANGCHECK);
1914 Check.FixTangent(Standard_True,Standard_True);
1916 IntTools_Curve aCurve;
1917 aCurve.SetCurve(BS);
1920 Handle(Geom2d_BSplineCurve) BS1=new Geom2d_BSplineCurve(tpoles2d,
1922 mbspc.Multiplicities(),
1924 GeomLib_Check2dBSplineCurve Check1(BS1,TOLCHECK,TOLANGCHECK);
1925 Check1.FixTangent(Standard_True,Standard_True);
1927 // ###########################################
1928 if(!rejectSurface && !reApprox) {
1929 Standard_Boolean isValid = IsCurveValid(BS1);
1931 reApprox = Standard_True;
1935 // ###########################################
1936 bPCurvesOk = CheckPCurve(BS1, myFace2);
1937 aCurve.SetSecondCurve2d(BS1);
1940 Handle(Geom2d_BSplineCurve) H2;
1941 aCurve.SetSecondCurve2d(H2);
1945 mbspc.Curve(1,tpoles2d);
1946 Handle(Geom2d_BSplineCurve) BS2=new Geom2d_BSplineCurve(tpoles2d,
1948 mbspc.Multiplicities(),
1950 GeomLib_Check2dBSplineCurve Check2(BS2,TOLCHECK,TOLANGCHECK);
1951 Check2.FixTangent(Standard_True,Standard_True);
1953 // ###########################################
1954 if(!rejectSurface && !reApprox) {
1955 Standard_Boolean isValid = IsCurveValid(BS2);
1957 reApprox = Standard_True;
1961 // ###########################################
1962 bPCurvesOk = bPCurvesOk && CheckPCurve(BS2, myFace1);
1963 aCurve.SetFirstCurve2d(BS2);
1966 Handle(Geom2d_BSplineCurve) H1;
1968 aCurve.SetFirstCurve2d(H1);
1971 //if points of the pcurves are out of the faces bounds
1972 //create 3d and 2d curves without approximation
1974 Handle(Geom2d_BSplineCurve) H1, H2;
1975 bPCurvesOk = Standard_True;
1977 Handle(Geom_Curve) aBSp=MakeBSpline(WL,ifprm, ilprm);
1980 H1 = MakeBSpline2d(WL, ifprm, ilprm, Standard_True);
1981 bPCurvesOk = CheckPCurve(H1, myFace1);
1985 H2 = MakeBSpline2d(WL, ifprm, ilprm, Standard_False);
1986 bPCurvesOk = bPCurvesOk && CheckPCurve(H2, myFace2);
1989 //if pcurves created without approximation are out of the
1990 //faces bounds, use approximated 3d and 2d curves
1992 IntTools_Curve aIC(aBSp, H1, H2);
1993 mySeqOfCurve.Append(aIC);
1995 mySeqOfCurve.Append(aCurve);
1998 mySeqOfCurve.Append(aCurve);
2001 }// else if(typs2 == GeomAbs_Plane)
2003 else { //typs2 != GeomAbs_Plane && typs1 != GeomAbs_Plane
2004 Standard_Boolean bIsValid1, bIsValid2;
2005 Handle(Geom_BSplineCurve) BS;
2006 Handle(Geom2d_BSplineCurve) aH2D;
2007 IntTools_Curve aCurve;
2009 bIsValid1=Standard_True;
2010 bIsValid2=Standard_True;
2012 const AppParCurves_MultiBSpCurve& mbspc = theapp3d.Value(j);
2013 nbpoles = mbspc.NbPoles();
2014 TColgp_Array1OfPnt tpoles(1,nbpoles);
2015 mbspc.Curve(1,tpoles);
2016 BS=new Geom_BSplineCurve(tpoles,
2018 mbspc.Multiplicities(),
2020 GeomLib_CheckBSplineCurve Check(BS,TOLCHECK,TOLANGCHECK);
2021 Check.FixTangent(Standard_True,Standard_True);
2023 aCurve.SetCurve(BS);
2024 aCurve.SetFirstCurve2d(aH2D);
2025 aCurve.SetSecondCurve2d(aH2D);
2029 Handle(Geom2d_BSplineCurve) BS1;
2030 TColgp_Array1OfPnt2d tpoles2d(1,nbpoles);
2031 mbspc.Curve(2,tpoles2d);
2033 BS1=new Geom2d_BSplineCurve(tpoles2d,
2035 mbspc.Multiplicities(),
2037 GeomLib_Check2dBSplineCurve newCheck(BS1,TOLCHECK,TOLANGCHECK);
2038 newCheck.FixTangent(Standard_True,Standard_True);
2041 bIsValid1=CheckPCurve(BS1, myFace1);
2044 aCurve.SetFirstCurve2d(BS1);
2047 Handle(Geom2d_BSplineCurve) BS1;
2048 fprm = BS->FirstParameter();
2049 lprm = BS->LastParameter();
2051 Handle(Geom2d_Curve) C2d;
2052 Standard_Real aTol = myTolApprox;
2053 GeomInt_IntSS::BuildPCurves(fprm, lprm, aTol,
2054 myHS1->ChangeSurface().Surface(), BS, C2d);
2055 BS1 = Handle(Geom2d_BSplineCurve)::DownCast(C2d);
2056 aCurve.SetFirstCurve2d(BS1);
2058 } // if(myApprox1) {
2062 Handle(Geom2d_BSplineCurve) BS2;
2063 TColgp_Array1OfPnt2d tpoles2d(1,nbpoles);
2064 mbspc.Curve((myApprox1==Standard_True)? 3 : 2,tpoles2d);
2065 BS2=new Geom2d_BSplineCurve(tpoles2d,
2067 mbspc.Multiplicities(),
2069 GeomLib_Check2dBSplineCurve newCheck(BS2,TOLCHECK,TOLANGCHECK);
2070 newCheck.FixTangent(Standard_True,Standard_True);
2073 bIsValid2=CheckPCurve(BS2, myFace2);
2075 aCurve.SetSecondCurve2d(BS2);
2078 Handle(Geom2d_BSplineCurve) BS2;
2079 fprm = BS->FirstParameter();
2080 lprm = BS->LastParameter();
2082 Handle(Geom2d_Curve) C2d;
2083 Standard_Real aTol = myTolApprox;
2084 GeomInt_IntSS::BuildPCurves(fprm, lprm, aTol,
2085 myHS2->ChangeSurface().Surface(), BS, C2d);
2086 BS2 = Handle(Geom2d_BSplineCurve)::DownCast(C2d);
2087 aCurve.SetSecondCurve2d(BS2);
2090 if (!bIsValid1 || !bIsValid2) {
2091 myTolApprox=aTolApproxImp;//1.e-5;
2092 tol2d = myTolApprox;
2093 reApprox = Standard_True;
2097 mySeqOfCurve.Append(aCurve);
2103 }// case IntPatch_Walking:{
2106 case IntPatch_Restriction:
2108 Handle(IntPatch_RLine) RL =
2109 Handle(IntPatch_RLine)::DownCast(L);
2110 Handle(Geom_Curve) aC3d;
2111 Handle(Geom2d_Curve) aC2d1, aC2d2;
2112 Standard_Real aTolReached;
2113 GeomInt_IntSS::TreatRLine(RL, myHS1, myHS2, aC3d,
2114 aC2d1, aC2d2, aTolReached);
2119 Bnd_Box2d aBox1, aBox2;
2121 const Standard_Real aU1f = myHS1->FirstUParameter(),
2122 aV1f = myHS1->FirstVParameter(),
2123 aU1l = myHS1->LastUParameter(),
2124 aV1l = myHS1->LastVParameter();
2125 const Standard_Real aU2f = myHS2->FirstUParameter(),
2126 aV2f = myHS2->FirstVParameter(),
2127 aU2l = myHS2->LastUParameter(),
2128 aV2l = myHS2->LastVParameter();
2130 aBox1.Add(gp_Pnt2d(aU1f, aV1f));
2131 aBox1.Add(gp_Pnt2d(aU1l, aV1l));
2132 aBox2.Add(gp_Pnt2d(aU2f, aV2f));
2133 aBox2.Add(gp_Pnt2d(aU2l, aV2l));
2135 GeomInt_VectorOfReal anArrayOfParameters;
2137 //We consider here that the intersection line is same-parameter-line
2138 anArrayOfParameters.Append(aC3d->FirstParameter());
2139 anArrayOfParameters.Append(aC3d->LastParameter());
2142 TrimILineOnSurfBoundaries(aC2d1, aC2d2, aBox1, aBox2, anArrayOfParameters);
2144 const Standard_Integer aNbIntersSolutionsm1 = anArrayOfParameters.Length() - 1;
2147 for(Standard_Integer anInd = 0; anInd < aNbIntersSolutionsm1; anInd++)
2149 const Standard_Real aParF = anArrayOfParameters(anInd),
2150 aParL = anArrayOfParameters(anInd+1);
2152 if((aParL - aParF) <= Precision::PConfusion())
2155 const Standard_Real aPar = 0.5*(aParF + aParL);
2158 Handle(Geom2d_Curve) aCurv2d1, aCurv2d2;
2161 aC2d1->D0(aPar, aPt);
2163 if(aBox1.IsOut(aPt))
2167 aCurv2d1 = new Geom2d_TrimmedCurve(aC2d1, aParF, aParL);
2172 aC2d2->D0(aPar, aPt);
2174 if(aBox2.IsOut(aPt))
2178 aCurv2d2 = new Geom2d_TrimmedCurve(aC2d2, aParF, aParL);
2181 Handle(Geom_Curve) aCurv3d = new Geom_TrimmedCurve(aC3d, aParF, aParL);
2183 IntTools_Curve aIC(aCurv3d, aCurv2d1, aCurv2d2);
2184 mySeqOfCurve.Append(aIC);
2194 //=======================================================================
2195 //function : Parameters
2197 //=======================================================================
2198 void Parameters(const Handle(GeomAdaptor_HSurface)& HS1,
2199 const Handle(GeomAdaptor_HSurface)& HS2,
2200 const gp_Pnt& Ptref,
2207 IntSurf_Quadric quad1,quad2;
2208 GeomAbs_SurfaceType typs = HS1->Surface().GetType();
2212 quad1.SetValue(HS1->Surface().Plane());
2214 case GeomAbs_Cylinder:
2215 quad1.SetValue(HS1->Surface().Cylinder());
2218 quad1.SetValue(HS1->Surface().Cone());
2220 case GeomAbs_Sphere:
2221 quad1.SetValue(HS1->Surface().Sphere());
2224 quad1.SetValue(HS1->Surface().Torus());
2227 Standard_ConstructionError::Raise("GeomInt_IntSS::MakeCurve");
2230 typs = HS2->Surface().GetType();
2233 quad2.SetValue(HS2->Surface().Plane());
2235 case GeomAbs_Cylinder:
2236 quad2.SetValue(HS2->Surface().Cylinder());
2239 quad2.SetValue(HS2->Surface().Cone());
2241 case GeomAbs_Sphere:
2242 quad2.SetValue(HS2->Surface().Sphere());
2245 quad2.SetValue(HS2->Surface().Torus());
2248 Standard_ConstructionError::Raise("GeomInt_IntSS::MakeCurve");
2251 quad1.Parameters(Ptref,U1,V1);
2252 quad2.Parameters(Ptref,U2,V2);
2255 //=======================================================================
2256 //function : MakeBSpline
2258 //=======================================================================
2259 Handle(Geom_Curve) MakeBSpline (const Handle(IntPatch_WLine)& WL,
2260 const Standard_Integer ideb,
2261 const Standard_Integer ifin)
2263 Standard_Integer i,nbpnt = ifin-ideb+1;
2264 TColgp_Array1OfPnt poles(1,nbpnt);
2265 TColStd_Array1OfReal knots(1,nbpnt);
2266 TColStd_Array1OfInteger mults(1,nbpnt);
2267 Standard_Integer ipidebm1;
2268 for(i=1,ipidebm1=i+ideb-1; i<=nbpnt;ipidebm1++, i++) {
2269 poles(i) = WL->Point(ipidebm1).Value();
2273 mults(1) = mults(nbpnt) = 2;
2275 new Geom_BSplineCurve(poles,knots,mults,1);
2279 //=======================================================================
2280 //function : MakeBSpline2d
2282 //=======================================================================
2283 Handle(Geom2d_BSplineCurve) MakeBSpline2d(const Handle(IntPatch_WLine)& theWLine,
2284 const Standard_Integer ideb,
2285 const Standard_Integer ifin,
2286 const Standard_Boolean onFirst)
2288 Standard_Integer i, nbpnt = ifin-ideb+1;
2289 TColgp_Array1OfPnt2d poles(1,nbpnt);
2290 TColStd_Array1OfReal knots(1,nbpnt);
2291 TColStd_Array1OfInteger mults(1,nbpnt);
2292 Standard_Integer ipidebm1;
2294 for(i = 1, ipidebm1 = i+ideb-1; i <= nbpnt; ipidebm1++, i++) {
2297 theWLine->Point(ipidebm1).ParametersOnS1(U, V);
2299 theWLine->Point(ipidebm1).ParametersOnS2(U, V);
2300 poles(i).SetCoord(U, V);
2304 mults(1) = mults(nbpnt) = 2;
2306 return new Geom2d_BSplineCurve(poles,knots,mults,1);
2308 //=======================================================================
2309 //function : PrepareLines3D
2311 //=======================================================================
2312 void IntTools_FaceFace::PrepareLines3D(const Standard_Boolean bToSplit)
2314 Standard_Integer i, aNbCurves;
2315 GeomAbs_SurfaceType aType1, aType2;
2316 IntTools_SequenceOfCurves aNewCvs;
2318 // 1. Treatment closed curves
2319 aNbCurves=mySeqOfCurve.Length();
2320 for (i=1; i<=aNbCurves; ++i) {
2321 const IntTools_Curve& aIC=mySeqOfCurve(i);
2324 Standard_Integer j, aNbC;
2325 IntTools_SequenceOfCurves aSeqCvs;
2327 aNbC=IntTools_Tools::SplitCurve(aIC, aSeqCvs);
2329 for (j=1; j<=aNbC; ++j) {
2330 const IntTools_Curve& aICNew=aSeqCvs(j);
2331 aNewCvs.Append(aICNew);
2335 aNewCvs.Append(aIC);
2339 aNewCvs.Append(aIC);
2343 // 2. Plane\Cone intersection when we had 4 curves
2344 aType1=myHS1->GetType();
2345 aType2=myHS2->GetType();
2346 aNbCurves=aNewCvs.Length();
2348 if ((aType1==GeomAbs_Plane && aType2==GeomAbs_Cone) ||
2349 (aType2==GeomAbs_Plane && aType1==GeomAbs_Cone)) {
2351 GeomAbs_CurveType aCType1;
2353 aCType1=aNewCvs(1).Type();
2354 if (aCType1==GeomAbs_Line) {
2355 IntTools_SequenceOfCurves aSeqIn, aSeqOut;
2357 for (i=1; i<=aNbCurves; ++i) {
2358 const IntTools_Curve& aIC=aNewCvs(i);
2362 IntTools_Tools::RejectLines(aSeqIn, aSeqOut);
2365 aNbCurves=aSeqOut.Length();
2366 for (i=1; i<=aNbCurves; ++i) {
2367 const IntTools_Curve& aIC=aSeqOut(i);
2368 aNewCvs.Append(aIC);
2372 }// if ((aType1==GeomAbs_Plane && aType2==GeomAbs_Cone)...
2374 // 3. Fill mySeqOfCurve
2375 mySeqOfCurve.Clear();
2376 aNbCurves=aNewCvs.Length();
2377 for (i=1; i<=aNbCurves; ++i) {
2378 const IntTools_Curve& aIC=aNewCvs(i);
2379 mySeqOfCurve.Append(aIC);
2382 //=======================================================================
2383 //function : CorrectSurfaceBoundaries
2385 //=======================================================================
2386 void CorrectSurfaceBoundaries(const TopoDS_Face& theFace,
2387 const Standard_Real theTolerance,
2388 Standard_Real& theumin,
2389 Standard_Real& theumax,
2390 Standard_Real& thevmin,
2391 Standard_Real& thevmax)
2393 Standard_Boolean enlarge, isuperiodic, isvperiodic;
2394 Standard_Real uinf, usup, vinf, vsup, delta;
2395 GeomAbs_SurfaceType aType;
2396 Handle(Geom_Surface) aSurface;
2398 aSurface = BRep_Tool::Surface(theFace);
2399 aSurface->Bounds(uinf, usup, vinf, vsup);
2400 delta = theTolerance;
2401 enlarge = Standard_False;
2403 GeomAdaptor_Surface anAdaptorSurface(aSurface);
2405 if(aSurface->IsKind(STANDARD_TYPE(Geom_RectangularTrimmedSurface))) {
2406 Handle(Geom_Surface) aBasisSurface =
2407 (Handle(Geom_RectangularTrimmedSurface)::DownCast(aSurface))->BasisSurface();
2409 if(aBasisSurface->IsKind(STANDARD_TYPE(Geom_RectangularTrimmedSurface)) ||
2410 aBasisSurface->IsKind(STANDARD_TYPE(Geom_OffsetSurface))) {
2415 if(aSurface->IsKind(STANDARD_TYPE(Geom_OffsetSurface))) {
2416 Handle(Geom_Surface) aBasisSurface =
2417 (Handle(Geom_OffsetSurface)::DownCast(aSurface))->BasisSurface();
2419 if(aBasisSurface->IsKind(STANDARD_TYPE(Geom_RectangularTrimmedSurface)) ||
2420 aBasisSurface->IsKind(STANDARD_TYPE(Geom_OffsetSurface))) {
2425 isuperiodic = anAdaptorSurface.IsUPeriodic();
2426 isvperiodic = anAdaptorSurface.IsVPeriodic();
2428 aType=anAdaptorSurface.GetType();
2429 if((aType==GeomAbs_BezierSurface) ||
2430 (aType==GeomAbs_BSplineSurface) ||
2431 (aType==GeomAbs_SurfaceOfExtrusion) ||
2432 (aType==GeomAbs_SurfaceOfRevolution) ||
2433 (aType==GeomAbs_Cylinder)) {
2434 enlarge=Standard_True;
2437 if(!isuperiodic && enlarge) {
2439 if(!Precision::IsInfinite(theumin) &&
2440 ((theumin - uinf) > delta))
2446 if(!Precision::IsInfinite(theumax) &&
2447 ((usup - theumax) > delta))
2453 if(!isvperiodic && enlarge) {
2454 if(!Precision::IsInfinite(thevmin) &&
2455 ((thevmin - vinf) > delta)) {
2461 if(!Precision::IsInfinite(thevmax) &&
2462 ((vsup - thevmax) > delta)) {
2471 Standard_Integer aNbP;
2472 Standard_Real aXP, dXfact, aXmid, aX1, aX2, aTolPA;
2474 aTolPA=Precision::Angular();
2477 aXP=anAdaptorSurface.UPeriod();
2478 dXfact=theumax-theumin;
2479 if (dXfact-aTolPA>aXP) {
2480 aXmid=0.5*(theumax+theumin);
2481 aNbP=RealToInt(aXmid/aXP);
2486 if (theumin>aTolPA) {
2487 aX1=theumin+aNbP*aXP;
2500 aXP=anAdaptorSurface.VPeriod();
2501 dXfact=thevmax-thevmin;
2502 if (dXfact-aTolPA>aXP) {
2503 aXmid=0.5*(thevmax+thevmin);
2504 aNbP=RealToInt(aXmid/aXP);
2509 if (thevmin>aTolPA) {
2510 aX1=thevmin+aNbP*aXP;
2523 if(isuperiodic || isvperiodic) {
2524 Standard_Boolean correct = Standard_False;
2525 Standard_Boolean correctU = Standard_False;
2526 Standard_Boolean correctV = Standard_False;
2528 TopExp_Explorer anExp;
2530 for(anExp.Init(theFace, TopAbs_EDGE); anExp.More(); anExp.Next()) {
2531 if(BRep_Tool::IsClosed(TopoDS::Edge(anExp.Current()), theFace)) {
2532 correct = Standard_True;
2534 TopoDS_Edge anEdge = TopoDS::Edge(anExp.Current());
2536 for(Standard_Integer i = 0; i < 2; i++) {
2538 anEdge.Orientation(TopAbs_FORWARD);
2541 anEdge.Orientation(TopAbs_REVERSED);
2543 Handle(Geom2d_Curve) aCurve = BRep_Tool::CurveOnSurface(anEdge, theFace, f, l);
2545 if(aCurve.IsNull()) {
2546 correct = Standard_False;
2549 Handle(Geom2d_Line) aLine = Handle(Geom2d_Line)::DownCast(aCurve);
2551 if(aLine.IsNull()) {
2552 correct = Standard_False;
2555 gp_Dir2d anUDir(1., 0.);
2556 gp_Dir2d aVDir(0., 1.);
2557 Standard_Real anAngularTolerance = Precision::Angular();
2559 correctU = correctU || aLine->Position().Direction().IsParallel(aVDir, anAngularTolerance);
2560 correctV = correctV || aLine->Position().Direction().IsParallel(anUDir, anAngularTolerance);
2562 gp_Pnt2d pp1 = aCurve->Value(f);
2564 gp_Pnt2d pp2 = aCurve->Value(l);
2573 Standard_Real umin, vmin, umax, vmax;
2574 aBox.Get(umin, vmin, umax, vmax);
2576 if(isuperiodic && correctU) {
2581 if(theumax > umax) {
2585 if(isvperiodic && correctV) {
2597 // The block is dedicated to determine whether WLine [ifprm, ilprm]
2598 // crosses the degenerated zone on each given surface or not.
2599 // If Yes -> We will not use info about surfaces during approximation
2600 // because inside degenerated zone of the surface the approx. algo.
2601 // uses wrong values of normal, etc., and resulting curve will have
2602 // oscillations that we would not like to have.
2607 Standard_Boolean IsDegeneratedZone(const gp_Pnt2d& aP2d,
2608 const Handle(Geom_Surface)& aS,
2609 const Standard_Integer iDir);
2611 Standard_Boolean IsPointInDegeneratedZone(const IntSurf_PntOn2S& aP2S,
2612 const TopoDS_Face& aF1,
2613 const TopoDS_Face& aF2);
2614 //=======================================================================
2615 //function : NotUseSurfacesForApprox
2617 //=======================================================================
2618 Standard_Boolean NotUseSurfacesForApprox(const TopoDS_Face& aF1,
2619 const TopoDS_Face& aF2,
2620 const Handle(IntPatch_WLine)& WL,
2621 const Standard_Integer ifprm,
2622 const Standard_Integer ilprm)
2624 Standard_Boolean bPInDZ;
2626 Handle(IntSurf_LineOn2S) aLineOn2S=WL->Curve();
2628 const IntSurf_PntOn2S& aP2Sfprm=aLineOn2S->Value(ifprm);
2629 bPInDZ=IsPointInDegeneratedZone(aP2Sfprm, aF1, aF2);
2634 const IntSurf_PntOn2S& aP2Slprm=aLineOn2S->Value(ilprm);
2635 bPInDZ=IsPointInDegeneratedZone(aP2Slprm, aF1, aF2);
2639 //=======================================================================
2640 //function : IsPointInDegeneratedZone
2642 //=======================================================================
2643 Standard_Boolean IsPointInDegeneratedZone(const IntSurf_PntOn2S& aP2S,
2644 const TopoDS_Face& aF1,
2645 const TopoDS_Face& aF2)
2648 Standard_Boolean bFlag=Standard_True;
2649 Standard_Real US11, US12, VS11, VS12, US21, US22, VS21, VS22;
2650 Standard_Real U1, V1, U2, V2, aDelta, aD;
2653 Handle(Geom_Surface)aS1 = BRep_Tool::Surface(aF1);
2654 aS1->Bounds(US11, US12, VS11, VS12);
2655 GeomAdaptor_Surface aGAS1(aS1);
2657 Handle(Geom_Surface)aS2 = BRep_Tool::Surface(aF2);
2658 aS1->Bounds(US21, US22, VS21, VS22);
2659 GeomAdaptor_Surface aGAS2(aS2);
2661 //const gp_Pnt& aP=aP2S.Value();
2662 aP2S.Parameters(U1, V1, U2, V2);
2666 aD=aGAS1.UResolution(aDelta);
2667 aP2d.SetCoord(U1, V1);
2668 if (fabs(U1-US11) < aD) {
2669 bFlag=IsDegeneratedZone(aP2d, aS1, 1);
2674 if (fabs(U1-US12) < aD) {
2675 bFlag=IsDegeneratedZone(aP2d, aS1, 1);
2680 aD=aGAS1.VResolution(aDelta);
2681 if (fabs(V1-VS11) < aDelta) {
2682 bFlag=IsDegeneratedZone(aP2d, aS1, 2);
2687 if (fabs(V1-VS12) < aDelta) {
2688 bFlag=IsDegeneratedZone(aP2d, aS1, 2);
2694 aD=aGAS2.UResolution(aDelta);
2695 aP2d.SetCoord(U2, V2);
2696 if (fabs(U2-US21) < aDelta) {
2697 bFlag=IsDegeneratedZone(aP2d, aS2, 1);
2702 if (fabs(U2-US22) < aDelta) {
2703 bFlag=IsDegeneratedZone(aP2d, aS2, 1);
2708 aD=aGAS2.VResolution(aDelta);
2709 if (fabs(V2-VS21) < aDelta) {
2710 bFlag=IsDegeneratedZone(aP2d, aS2, 2);
2715 if (fabs(V2-VS22) < aDelta) {
2716 bFlag=IsDegeneratedZone(aP2d, aS2, 2);
2724 //=======================================================================
2725 //function : IsDegeneratedZone
2727 //=======================================================================
2728 Standard_Boolean IsDegeneratedZone(const gp_Pnt2d& aP2d,
2729 const Handle(Geom_Surface)& aS,
2730 const Standard_Integer iDir)
2732 Standard_Boolean bFlag=Standard_True;
2733 Standard_Real US1, US2, VS1, VS2, dY, dX, d1, d2, dD;
2734 Standard_Real aXm, aYm, aXb, aYb, aXe, aYe;
2735 aS->Bounds(US1, US2, VS1, VS2);
2737 gp_Pnt aPm, aPb, aPe;
2742 aS->D0(aXm, aYm, aPm);
2759 aS->D0(aXb, aYb, aPb);
2760 aS->D0(aXe, aYe, aPe);
2762 d1=aPm.Distance(aPb);
2763 d2=aPm.Distance(aPe);
2764 if (d1 < dD && d2 < dD) {
2781 aS->D0(aXb, aYb, aPb);
2782 aS->D0(aXe, aYe, aPe);
2784 d1=aPm.Distance(aPb);
2785 d2=aPm.Distance(aPe);
2786 if (d1 < dD && d2 < dD) {
2794 // Check if aNextPnt lies inside of tube build on aBasePnt and aBaseVec.
2796 static Standard_Boolean IsInsideIn2d(const gp_Pnt2d& aBasePnt,
2797 const gp_Vec2d& aBaseVec,
2798 const gp_Pnt2d& aNextPnt,
2799 const Standard_Real aSquareMaxDist)
2801 gp_Vec2d aVec2d(aBasePnt, aNextPnt);
2803 //d*d = (basevec^(nextpnt-basepnt))**2 / basevec**2
2804 Standard_Real aCross = aVec2d.Crossed(aBaseVec);
2805 Standard_Real aSquareDist = aCross * aCross
2806 / aBaseVec.SquareMagnitude();
2808 return (aSquareDist <= aSquareMaxDist);
2811 // Check if aNextPnt lies inside of tube build on aBasePnt and aBaseVec.
2813 static Standard_Boolean IsInsideIn3d(const gp_Pnt& aBasePnt,
2814 const gp_Vec& aBaseVec,
2815 const gp_Pnt& aNextPnt,
2816 const Standard_Real aSquareMaxDist)
2818 gp_Vec aVec(aBasePnt, aNextPnt);
2820 //d*d = (basevec^(nextpnt-basepnt))**2 / basevec**2
2821 Standard_Real aSquareDist = aVec.CrossSquareMagnitude(aBaseVec)
2822 / aBaseVec.SquareMagnitude();
2824 return (aSquareDist <= aSquareMaxDist);
2827 //=========================================================================
2828 // static function : ComputePurgedWLine
2829 // purpose : Removes equal points (leave one of equal points) from theWLine
2830 // and recompute vertex parameters.
2831 // Removes exceed points using tube criteria:
2832 // delete 7D point if it lies near to expected lines in 2d and 3d.
2833 // Each task (2d, 2d, 3d) have its own tolerance and checked separately.
2834 // Returns new WLine or null WLine if the number
2835 // of the points is less than 2.
2836 //=========================================================================
2837 Handle(IntPatch_WLine) ComputePurgedWLine(const Handle(IntPatch_WLine) &theWLine,
2838 const Handle(GeomAdaptor_HSurface) &theS1,
2839 const Handle(GeomAdaptor_HSurface) &theS2)
2841 Standard_Integer i, k, v, nb, nbvtx;
2842 Handle(IntPatch_WLine) aResult;
2843 nbvtx = theWLine->NbVertex();
2844 nb = theWLine->NbPnts();
2846 const IntSurf_PntOn2S& p1 = theWLine->Point(1);
2847 const IntSurf_PntOn2S& p2 = theWLine->Point(2);
2848 if(p1.Value().IsEqual(p2.Value(), gp::Resolution())) {
2853 Handle(IntPatch_WLine) aLocalWLine;
2854 Handle(IntPatch_WLine) aTmpWLine = theWLine;
2855 Handle(IntSurf_LineOn2S) aLineOn2S = new IntSurf_LineOn2S();
2856 aLocalWLine = new IntPatch_WLine(aLineOn2S, Standard_False);
2857 for(i = 1; i <= nb; i++) {
2858 aLineOn2S->Add(theWLine->Point(i));
2861 for(v = 1; v <= nbvtx; v++) {
2862 aLocalWLine->AddVertex(theWLine->Vertex(v));
2864 for(i = 1; i <= aLineOn2S->NbPoints(); i++) {
2865 Standard_Integer aStartIndex = i + 1;
2866 Standard_Integer anEndIndex = i + 5;
2867 nb = aLineOn2S->NbPoints();
2868 anEndIndex = (anEndIndex > nb) ? nb : anEndIndex;
2870 if((aStartIndex > nb) || (anEndIndex <= 1)) {
2875 while(k <= anEndIndex) {
2878 IntSurf_PntOn2S p1 = aLineOn2S->Value(i);
2879 IntSurf_PntOn2S p2 = aLineOn2S->Value(k);
2881 Standard_Real UV[8];
2882 p1.Parameters(UV[0], UV[1], UV[2], UV[3]);
2883 p2.Parameters(UV[4], UV[5], UV[6], UV[7]);
2885 Standard_Real aMax = Abs(UV[0]);
2886 for(Standard_Integer anIdx = 1; anIdx < 8; anIdx++)
2888 if (aMax < Abs(UV[anIdx]))
2889 aMax = Abs(UV[anIdx]);
2892 if(p1.Value().IsEqual(p2.Value(), gp::Resolution()) ||
2893 Abs(UV[0] - UV[4]) + Abs(UV[1] - UV[5]) < 1.0e-16 * aMax ||
2894 Abs(UV[2] - UV[6]) + Abs(UV[3] - UV[7]) < 1.0e-16 * aMax )
2896 aTmpWLine = aLocalWLine;
2897 aLocalWLine = new IntPatch_WLine(aLineOn2S, Standard_False);
2899 for(v = 1; v <= aTmpWLine->NbVertex(); v++) {
2900 IntPatch_Point aVertex = aTmpWLine->Vertex(v);
2901 Standard_Integer avertexindex = (Standard_Integer)aVertex.ParameterOnLine();
2903 if(avertexindex >= k) {
2904 aVertex.SetParameter(aVertex.ParameterOnLine() - 1.);
2906 aLocalWLine->AddVertex(aVertex);
2908 aLineOn2S->RemovePoint(k);
2917 if (aLineOn2S->NbPoints() <= 2)
2919 if (aLineOn2S->NbPoints() == 2)
2925 const Standard_Integer aMinNbBadDistr = 15;
2926 const Standard_Integer aNbSingleBezier = 30;
2927 // Avoid purge in case of C0 continuity:
2928 // Intersection approximator may produce invalid curve after purge, example:
2929 // bugs modalg_5 bug24731,
2930 // Do not run purger when base number of points is too small.
2931 if (theS1->UContinuity() == GeomAbs_C0 ||
2932 theS1->VContinuity() == GeomAbs_C0 ||
2933 theS2->UContinuity() == GeomAbs_C0 ||
2934 theS2->VContinuity() == GeomAbs_C0 ||
2935 nb < aNbSingleBezier)
2940 // 1 - Delete point.
2942 // -1 - Vertex point (not delete).
2943 NCollection_Array1<Standard_Integer> aNewPointsHash(1, aLineOn2S->NbPoints());
2944 for(i = 1; i <= aLineOn2S->NbPoints(); i++)
2945 aNewPointsHash.SetValue(i, 0);
2947 for(v = 1; v <= aLocalWLine->NbVertex(); v++)
2949 IntPatch_Point aVertex = aLocalWLine->Vertex(v);
2950 Standard_Integer avertexindex = (Standard_Integer)aVertex.ParameterOnLine();
2951 aNewPointsHash.SetValue(avertexindex, -1);
2954 // Workaround to handle case of small amount points after purge.
2955 // Test "boolean boptuc_complex B5" and similar.
2956 Standard_Integer aNbPnt = 0;
2958 // Inital computations.
2959 Standard_Real UonS1[3], VonS1[3], UonS2[3], VonS2[3];
2960 aLineOn2S->Value(1).ParametersOnS1(UonS1[0], VonS1[0]);
2961 aLineOn2S->Value(2).ParametersOnS1(UonS1[1], VonS1[1]);
2962 aLineOn2S->Value(1).ParametersOnS2(UonS2[0], VonS2[0]);
2963 aLineOn2S->Value(2).ParametersOnS2(UonS2[1], VonS2[1]);
2965 gp_Pnt2d aBase2dPnt1(UonS1[0], VonS1[0]);
2966 gp_Pnt2d aBase2dPnt2(UonS2[0], VonS2[0]);
2967 gp_Vec2d aBase2dVec1(UonS1[1] - UonS1[0], VonS1[1] - VonS1[0]);
2968 gp_Vec2d aBase2dVec2(UonS2[1] - UonS2[0], VonS2[1] - VonS2[0]);
2969 gp_Pnt aBase3dPnt = aLineOn2S->Value(1).Value();
2970 gp_Vec aBase3dVec(aLineOn2S->Value(1).Value(), aLineOn2S->Value(2).Value());
2972 // Choose base tolerance and scale it to pipe algorithm.
2973 const Standard_Real aBaseTolerance = Precision::Approximation();
2974 Standard_Real aResS1Tol = Min(theS1->UResolution(aBaseTolerance),
2975 theS1->VResolution(aBaseTolerance));
2976 Standard_Real aResS2Tol = Min(theS2->UResolution(aBaseTolerance),
2977 theS2->VResolution(aBaseTolerance));
2978 Standard_Real aTol1 = aResS1Tol * aResS1Tol;
2979 Standard_Real aTol2 = aResS2Tol * aResS2Tol;
2980 Standard_Real aTol3d = aBaseTolerance * aBaseTolerance;
2982 const Standard_Real aLimitCoeff = 0.99 * 0.99;
2983 for(i = 3; i <= aLineOn2S->NbPoints(); i++)
2985 Standard_Boolean isDeleteState = Standard_False;
2987 aLineOn2S->Value(i).ParametersOnS1(UonS1[2], VonS1[2]);
2988 aLineOn2S->Value(i).ParametersOnS2(UonS2[2], VonS2[2]);
2989 gp_Pnt2d aPnt2dOnS1(UonS1[2], VonS1[2]);
2990 gp_Pnt2d aPnt2dOnS2(UonS2[2], VonS2[2]);
2991 const gp_Pnt& aPnt3d = aLineOn2S->Value(i).Value();
2993 if (aNewPointsHash(i - 1) != - 1 &&
2994 IsInsideIn2d(aBase2dPnt1, aBase2dVec1, aPnt2dOnS1, aTol1) &&
2995 IsInsideIn2d(aBase2dPnt2, aBase2dVec2, aPnt2dOnS2, aTol2) &&
2996 IsInsideIn3d(aBase3dPnt, aBase3dVec, aPnt3d, aTol3d) )
2998 // Handle possible uneven parametrization on one of 2d subspaces.
2999 // Delete point only when expected lengths are close to each other (aLimitCoeff).
3003 // c2d2 - geometrically line, but have uneven parametrization -> c2d2 is bspline.
3004 gp_XY aPntOnS1[2]= { gp_XY(UonS1[1] - UonS1[0], VonS1[1] - VonS1[0])
3005 , gp_XY(UonS1[2] - UonS1[1], VonS1[2] - VonS1[1])};
3006 gp_XY aPntOnS2[2]= { gp_XY(UonS2[1] - UonS2[0], VonS2[1] - VonS2[0])
3007 , gp_XY(UonS2[2] - UonS2[1], VonS2[2] - VonS2[1])};
3009 Standard_Real aStepOnS1 = aPntOnS1[0].SquareModulus() / aPntOnS1[1].SquareModulus();
3010 Standard_Real aStepOnS2 = aPntOnS2[0].SquareModulus() / aPntOnS2[1].SquareModulus();
3012 Standard_Real aStepCoeff = Min(aStepOnS1, aStepOnS2) / Max(aStepOnS1, aStepOnS2);
3014 if (aStepCoeff > aLimitCoeff)
3016 // Set hash flag to "Delete" state.
3017 isDeleteState = Standard_True;
3018 aNewPointsHash.SetValue(i - 1, 1);
3020 // Change middle point.
3021 UonS1[1] = UonS1[2];
3022 UonS2[1] = UonS2[2];
3023 VonS1[1] = VonS1[2];
3024 VonS2[1] = VonS2[2];
3030 // Compute new pipe parameters.
3031 UonS1[0] = UonS1[1];
3032 VonS1[0] = VonS1[1];
3033 UonS2[0] = UonS2[1];
3034 VonS2[0] = VonS2[1];
3036 UonS1[1] = UonS1[2];
3037 VonS1[1] = VonS1[2];
3038 UonS2[1] = UonS2[2];
3039 VonS2[1] = VonS2[2];
3041 aBase2dPnt1.SetCoord(UonS1[0], VonS1[0]);
3042 aBase2dPnt2.SetCoord(UonS2[0], VonS2[0]);
3043 aBase2dVec1.SetCoord(UonS1[1] - UonS1[0], VonS1[1] - VonS1[0]);
3044 aBase2dVec2.SetCoord(UonS2[1] - UonS2[0], VonS2[1] - VonS2[0]);
3045 aBase3dPnt = aLineOn2S->Value(i - 1).Value();
3046 aBase3dVec = gp_Vec(aLineOn2S->Value(i - 1).Value(), aLineOn2S->Value(i).Value());
3052 // Workaround to handle case of small amount of points after purge.
3053 // Test "boolean boptuc_complex B5" and similar.
3054 // This is possible since there are at least two points.
3055 if (aNewPointsHash(1) == -1 &&
3056 aNewPointsHash(2) == -1 &&
3060 aNewPointsHash(1) = 1;
3062 if (aNewPointsHash(aLineOn2S->NbPoints() - 1) == -1 &&
3063 aNewPointsHash(aLineOn2S->NbPoints() ) == -1 &&
3067 aNewPointsHash(aLineOn2S->NbPoints() ) = 1;
3070 // Purgre when too small amount of points left.
3073 for(i = aNewPointsHash.Lower(); i <= aNewPointsHash.Upper(); i++)
3075 if (aNewPointsHash(i) != -1)
3077 aNewPointsHash(i) = 1;
3082 // Handle possible bad distribution of points,
3083 // which are will converted into one single bezier curve (less than 30 points).
3084 // Make distribution more even:
3085 // max step will be nearly to 0.1 of param distance.
3086 if (aNbPnt + 2 > aMinNbBadDistr &&
3087 aNbPnt + 2 < aNbSingleBezier )
3089 for(Standard_Integer anIdx = 1; anIdx <= 8; anIdx++)
3091 Standard_Integer aHashIdx =
3092 Standard_Integer(anIdx * aLineOn2S->NbPoints() / 9);
3095 aNewPointsHash(aHashIdx) = 0;
3099 aTmpWLine = aLocalWLine;
3100 Handle(IntSurf_LineOn2S) aPurgedLineOn2S = new IntSurf_LineOn2S();
3101 aLocalWLine = new IntPatch_WLine(aPurgedLineOn2S, Standard_False);
3102 Standard_Integer anOldLineIdx = 1, aVertexIdx = 1;
3103 for(i = 1; i <= aNewPointsHash.Upper(); i++)
3105 if (aNewPointsHash(i) == 0)
3107 // Store this point.
3108 aPurgedLineOn2S->Add(aLineOn2S->Value(i));
3111 else if (aNewPointsHash(i) == -1)
3114 IntPatch_Point aVertex = aTmpWLine->Vertex(aVertexIdx++);
3115 aVertex.SetParameter(anOldLineIdx++);
3116 aLocalWLine->AddVertex(aVertex);
3117 aPurgedLineOn2S->Add(aLineOn2S->Value(i));
3121 if(aPurgedLineOn2S->NbPoints() > 1) {
3122 aResult = aLocalWLine;
3127 //=======================================================================
3130 //=======================================================================
3131 void TolR3d(const TopoDS_Face& aF1,
3132 const TopoDS_Face& aF2,
3133 Standard_Real& myTolReached3d)
3135 Standard_Real aTolF1, aTolF2, aTolFMax, aTolTresh;
3137 aTolTresh=2.999999e-3;
3138 aTolF1 = BRep_Tool::Tolerance(aF1);
3139 aTolF2 = BRep_Tool::Tolerance(aF2);
3140 aTolFMax=Max(aTolF1, aTolF2);
3142 if (aTolFMax>aTolTresh) {
3143 myTolReached3d=aTolFMax;
3146 //=======================================================================
3147 //function : IsPointOnBoundary
3149 //=======================================================================
3150 Standard_Boolean IsPointOnBoundary(const Standard_Real theParameter,
3151 const Standard_Real theFirstBoundary,
3152 const Standard_Real theSecondBoundary,
3153 const Standard_Real theResolution,
3154 Standard_Boolean& IsOnFirstBoundary)
3156 Standard_Boolean bRet;
3158 Standard_Real adist;
3160 bRet=Standard_False;
3161 for(i = 0; i < 2; ++i) {
3162 IsOnFirstBoundary = (i == 0);
3163 if (IsOnFirstBoundary) {
3164 adist = fabs(theParameter - theFirstBoundary);
3167 adist = fabs(theParameter - theSecondBoundary);
3169 if(adist < theResolution) {
3175 // ------------------------------------------------------------------------------------------------
3176 // static function: FindPoint
3178 // ------------------------------------------------------------------------------------------------
3179 Standard_Boolean FindPoint(const gp_Pnt2d& theFirstPoint,
3180 const gp_Pnt2d& theLastPoint,
3181 const Standard_Real theUmin,
3182 const Standard_Real theUmax,
3183 const Standard_Real theVmin,
3184 const Standard_Real theVmax,
3185 gp_Pnt2d& theNewPoint) {
3187 gp_Vec2d aVec(theFirstPoint, theLastPoint);
3188 Standard_Integer i = 0, j = 0;
3190 for(i = 0; i < 4; i++) {
3191 gp_Vec2d anOtherVec;
3192 gp_Vec2d anOtherVecNormal;
3193 gp_Pnt2d aprojpoint = theLastPoint;
3196 anOtherVec.SetX(0.);
3197 anOtherVec.SetY(1.);
3198 anOtherVecNormal.SetX(1.);
3199 anOtherVecNormal.SetY(0.);
3202 aprojpoint.SetX(theUmin);
3204 aprojpoint.SetX(theUmax);
3207 anOtherVec.SetX(1.);
3208 anOtherVec.SetY(0.);
3209 anOtherVecNormal.SetX(0.);
3210 anOtherVecNormal.SetY(1.);
3213 aprojpoint.SetY(theVmin);
3215 aprojpoint.SetY(theVmax);
3217 gp_Vec2d anormvec = aVec;
3218 anormvec.Normalize();
3219 RefineVector(anormvec);
3220 Standard_Real adot1 = anormvec.Dot(anOtherVecNormal);
3222 if(fabs(adot1) < Precision::Angular())
3224 Standard_Real adist = 0.;
3225 Standard_Boolean bIsOut = Standard_False;
3228 adist = (i < 2) ? fabs(theLastPoint.X() - theUmin) : fabs(theLastPoint.X() - theUmax);
3229 bIsOut = (i < 2) ? (theLastPoint.X() < theUmin) : (theLastPoint.X() > theUmax);
3232 adist = (i < 2) ? fabs(theLastPoint.Y() - theVmin) : fabs(theLastPoint.Y() - theVmax);
3233 bIsOut = (i < 2) ? (theLastPoint.Y() < theVmin) : (theLastPoint.Y() > theVmax);
3235 Standard_Real anoffset = adist * anOtherVec.Dot(anormvec) / adot1;
3237 for(j = 0; j < 2; j++) {
3238 anoffset = (j == 0) ? anoffset : -anoffset;
3239 gp_Pnt2d acurpoint(aprojpoint.XY() + (anOtherVec.XY()*anoffset));
3240 gp_Vec2d acurvec(theLastPoint, acurpoint);
3244 Standard_Real aDotX, anAngleX;
3246 aDotX = aVec.Dot(acurvec);
3247 anAngleX = aVec.Angle(acurvec);
3249 if(aDotX > 0. && fabs(anAngleX) < Precision::PConfusion()) {
3251 if((acurpoint.Y() >= theVmin) &&
3252 (acurpoint.Y() <= theVmax)) {
3253 theNewPoint = acurpoint;
3254 return Standard_True;
3258 if((acurpoint.X() >= theUmin) &&
3259 (acurpoint.X() <= theUmax)) {
3260 theNewPoint = acurpoint;
3261 return Standard_True;
3267 return Standard_False;
3271 // ------------------------------------------------------------------------------------------------
3272 // static function: FindPoint
3273 // purpose: Find point on the boundary of radial tangent zone
3274 // ------------------------------------------------------------------------------------------------
3275 Standard_Boolean FindPoint(const gp_Pnt2d& theFirstPoint,
3276 const gp_Pnt2d& theLastPoint,
3277 const Standard_Real theUmin,
3278 const Standard_Real theUmax,
3279 const Standard_Real theVmin,
3280 const Standard_Real theVmax,
3281 const gp_Pnt2d& theTanZoneCenter,
3282 const Standard_Real theZoneRadius,
3283 Handle(GeomAdaptor_HSurface) theGASurface,
3284 gp_Pnt2d& theNewPoint) {
3285 theNewPoint = theLastPoint;
3287 if ( !IsInsideTanZone( theLastPoint, theTanZoneCenter, theZoneRadius, theGASurface) )
3288 return Standard_False;
3290 Standard_Real aUResolution = theGASurface->UResolution( theZoneRadius );
3291 Standard_Real aVResolution = theGASurface->VResolution( theZoneRadius );
3293 Standard_Real aRadius = ( aUResolution < aVResolution ) ? aUResolution : aVResolution;
3294 gp_Ax22d anAxis( theTanZoneCenter, gp_Dir2d(1, 0), gp_Dir2d(0, 1) );
3295 gp_Circ2d aCircle( anAxis, aRadius );
3298 gp_Vec2d aDir( theLastPoint.XY() - theFirstPoint.XY() );
3299 Standard_Real aLength = aDir.Magnitude();
3300 if ( aLength <= gp::Resolution() )
3301 return Standard_False;
3302 gp_Lin2d aLine( theFirstPoint, aDir );
3305 Handle(Geom2d_Line) aCLine = new Geom2d_Line( aLine );
3306 Handle(Geom2d_TrimmedCurve) aC1 = new Geom2d_TrimmedCurve( aCLine, 0, aLength );
3307 Handle(Geom2d_Circle) aC2 = new Geom2d_Circle( aCircle );
3309 Standard_Real aTol = aRadius * 0.001;
3310 aTol = ( aTol < Precision::PConfusion() ) ? Precision::PConfusion() : aTol;
3312 Geom2dAPI_InterCurveCurve anIntersector;
3313 anIntersector.Init( aC1, aC2, aTol );
3315 if ( anIntersector.NbPoints() == 0 )
3316 return Standard_False;
3318 Standard_Boolean aFound = Standard_False;
3319 Standard_Real aMinDist = aLength * aLength;
3320 Standard_Integer i = 0;
3321 for ( i = 1; i <= anIntersector.NbPoints(); i++ ) {
3322 gp_Pnt2d aPInt = anIntersector.Point( i );
3323 if ( aPInt.SquareDistance( theFirstPoint ) < aMinDist ) {
3324 if ( ( aPInt.X() >= theUmin ) && ( aPInt.X() <= theUmax ) &&
3325 ( aPInt.Y() >= theVmin ) && ( aPInt.Y() <= theVmax ) ) {
3326 theNewPoint = aPInt;
3327 aFound = Standard_True;
3335 // ------------------------------------------------------------------------------------------------
3336 // static function: IsInsideTanZone
3337 // purpose: Check if point is inside a radial tangent zone
3338 // ------------------------------------------------------------------------------------------------
3339 Standard_Boolean IsInsideTanZone(const gp_Pnt2d& thePoint,
3340 const gp_Pnt2d& theTanZoneCenter,
3341 const Standard_Real theZoneRadius,
3342 Handle(GeomAdaptor_HSurface) theGASurface) {
3344 Standard_Real aUResolution = theGASurface->UResolution( theZoneRadius );
3345 Standard_Real aVResolution = theGASurface->VResolution( theZoneRadius );
3346 Standard_Real aRadiusSQR = ( aUResolution < aVResolution ) ? aUResolution : aVResolution;
3347 aRadiusSQR *= aRadiusSQR;
3348 if ( thePoint.SquareDistance( theTanZoneCenter ) <= aRadiusSQR )
3349 return Standard_True;
3350 return Standard_False;
3353 // ------------------------------------------------------------------------------------------------
3354 // static function: CheckTangentZonesExist
3355 // purpose: Check if tangent zone exists
3356 // ------------------------------------------------------------------------------------------------
3357 Standard_Boolean CheckTangentZonesExist( const Handle(GeomAdaptor_HSurface)& theSurface1,
3358 const Handle(GeomAdaptor_HSurface)& theSurface2 )
3360 if ( ( theSurface1->GetType() != GeomAbs_Torus ) ||
3361 ( theSurface2->GetType() != GeomAbs_Torus ) )
3362 return Standard_False;
3364 gp_Torus aTor1 = theSurface1->Torus();
3365 gp_Torus aTor2 = theSurface2->Torus();
3367 if ( aTor1.Location().Distance( aTor2.Location() ) > Precision::Confusion() )
3368 return Standard_False;
3370 if ( ( fabs( aTor1.MajorRadius() - aTor2.MajorRadius() ) > Precision::Confusion() ) ||
3371 ( fabs( aTor1.MinorRadius() - aTor2.MinorRadius() ) > Precision::Confusion() ) )
3372 return Standard_False;
3374 if ( ( aTor1.MajorRadius() < aTor1.MinorRadius() ) ||
3375 ( aTor2.MajorRadius() < aTor2.MinorRadius() ) )
3376 return Standard_False;
3377 return Standard_True;
3380 // ------------------------------------------------------------------------------------------------
3381 // static function: ComputeTangentZones
3383 // ------------------------------------------------------------------------------------------------
3384 Standard_Integer ComputeTangentZones( const Handle(GeomAdaptor_HSurface)& theSurface1,
3385 const Handle(GeomAdaptor_HSurface)& theSurface2,
3386 const TopoDS_Face& theFace1,
3387 const TopoDS_Face& theFace2,
3388 Handle(TColgp_HArray1OfPnt2d)& theResultOnS1,
3389 Handle(TColgp_HArray1OfPnt2d)& theResultOnS2,
3390 Handle(TColStd_HArray1OfReal)& theResultRadius,
3391 const Handle(IntTools_Context)& aContext)
3393 Standard_Integer aResult = 0;
3394 if ( !CheckTangentZonesExist( theSurface1, theSurface2 ) )
3398 TColgp_SequenceOfPnt2d aSeqResultS1, aSeqResultS2;
3399 TColStd_SequenceOfReal aSeqResultRad;
3401 gp_Torus aTor1 = theSurface1->Torus();
3402 gp_Torus aTor2 = theSurface2->Torus();
3404 gp_Ax2 anax1( aTor1.Location(), aTor1.Axis().Direction() );
3405 gp_Ax2 anax2( aTor2.Location(), aTor2.Axis().Direction() );
3406 Standard_Integer j = 0;
3408 for ( j = 0; j < 2; j++ ) {
3409 Standard_Real aCoef = ( j == 0 ) ? -1 : 1;
3410 Standard_Real aRadius1 = fabs(aTor1.MajorRadius() + aCoef * aTor1.MinorRadius());
3411 Standard_Real aRadius2 = fabs(aTor2.MajorRadius() + aCoef * aTor2.MinorRadius());
3413 gp_Circ aCircle1( anax1, aRadius1 );
3414 gp_Circ aCircle2( anax2, aRadius2 );
3416 // roughly compute radius of tangent zone for perpendicular case
3417 Standard_Real aCriteria = Precision::Confusion() * 0.5;
3419 Standard_Real aT1 = aCriteria;
3420 Standard_Real aT2 = aCriteria;
3422 // internal tangency
3423 Standard_Real aR = ( aRadius1 > aTor2.MinorRadius() ) ? aRadius1 : aTor2.MinorRadius();
3424 //aT1 = aCriteria * aCriteria + aR * aR - ( aR - aCriteria ) * ( aR - aCriteria );
3425 aT1 = 2. * aR * aCriteria;
3429 // external tangency
3430 Standard_Real aRb = ( aRadius1 > aTor2.MinorRadius() ) ? aRadius1 : aTor2.MinorRadius();
3431 Standard_Real aRm = ( aRadius1 < aTor2.MinorRadius() ) ? aRadius1 : aTor2.MinorRadius();
3432 Standard_Real aDelta = aRb - aCriteria;
3434 aDelta -= aRm * aRm;
3435 aDelta /= 2. * (aRb - aRm);
3436 aDelta -= 0.5 * (aRb - aRm);
3438 aT1 = 2. * aRm * (aRm - aDelta);
3441 aCriteria = ( aT1 > aT2) ? aT1 : aT2;
3442 if ( aCriteria > 0 )
3443 aCriteria = sqrt( aCriteria );
3445 if ( aCriteria > 0.5 * aTor1.MinorRadius() ) {
3446 // too big zone -> drop to minimum
3447 aCriteria = Precision::Confusion();
3450 GeomAdaptor_Curve aC1( new Geom_Circle(aCircle1) );
3451 GeomAdaptor_Curve aC2( new Geom_Circle(aCircle2) );
3452 Extrema_ExtCC anExtrema(aC1, aC2, 0, 2. * M_PI, 0, 2. * M_PI,
3453 Precision::PConfusion(), Precision::PConfusion());
3455 if ( anExtrema.IsDone() ) {
3457 Standard_Integer i = 0;
3458 for ( i = 1; i <= anExtrema.NbExt(); i++ ) {
3459 if ( anExtrema.SquareDistance(i) > aCriteria * aCriteria )
3462 Extrema_POnCurv P1, P2;
3463 anExtrema.Points( i, P1, P2 );
3465 Standard_Boolean bFoundResult = Standard_True;
3468 Standard_Integer surfit = 0;
3469 for ( surfit = 0; surfit < 2; surfit++ ) {
3470 GeomAPI_ProjectPointOnSurf& aProjector =
3471 (surfit == 0) ? aContext->ProjPS(theFace1) : aContext->ProjPS(theFace2);
3473 gp_Pnt aP3d = (surfit == 0) ? P1.Value() : P2.Value();
3474 aProjector.Perform(aP3d);
3476 if(!aProjector.IsDone())
3477 bFoundResult = Standard_False;
3479 if(aProjector.LowerDistance() > aCriteria) {
3480 bFoundResult = Standard_False;
3483 Standard_Real foundU = 0, foundV = 0;
3484 aProjector.LowerDistanceParameters(foundU, foundV);
3486 pr1 = gp_Pnt2d( foundU, foundV );
3488 pr2 = gp_Pnt2d( foundU, foundV );
3492 if ( bFoundResult ) {
3493 aSeqResultS1.Append( pr1 );
3494 aSeqResultS2.Append( pr2 );
3495 aSeqResultRad.Append( aCriteria );
3497 // torus is u and v periodic
3498 const Standard_Real twoPI = M_PI + M_PI;
3499 Standard_Real arr1tmp[2] = {pr1.X(), pr1.Y()};
3500 Standard_Real arr2tmp[2] = {pr2.X(), pr2.Y()};
3502 // iteration on period bounds
3503 for ( Standard_Integer k1 = 0; k1 < 2; k1++ ) {
3504 Standard_Real aBound = ( k1 == 0 ) ? 0 : twoPI;
3505 Standard_Real aShift = ( k1 == 0 ) ? twoPI : -twoPI;
3507 // iteration on surfaces
3508 for ( Standard_Integer k2 = 0; k2 < 2; k2++ ) {
3509 Standard_Real* arr1 = ( k2 == 0 ) ? arr1tmp : arr2tmp;
3510 Standard_Real* arr2 = ( k2 != 0 ) ? arr1tmp : arr2tmp;
3511 TColgp_SequenceOfPnt2d& aSeqS1 = ( k2 == 0 ) ? aSeqResultS1 : aSeqResultS2;
3512 TColgp_SequenceOfPnt2d& aSeqS2 = ( k2 != 0 ) ? aSeqResultS1 : aSeqResultS2;
3514 if (fabs(arr1[0] - aBound) < Precision::PConfusion()) {
3515 aSeqS1.Append( gp_Pnt2d( arr1[0] + aShift, arr1[1] ) );
3516 aSeqS2.Append( gp_Pnt2d( arr2[0], arr2[1] ) );
3517 aSeqResultRad.Append( aCriteria );
3519 if (fabs(arr1[1] - aBound) < Precision::PConfusion()) {
3520 aSeqS1.Append( gp_Pnt2d( arr1[0], arr1[1] + aShift) );
3521 aSeqS2.Append( gp_Pnt2d( arr2[0], arr2[1] ) );
3522 aSeqResultRad.Append( aCriteria );
3530 aResult = aSeqResultRad.Length();
3532 if ( aResult > 0 ) {
3533 theResultOnS1 = new TColgp_HArray1OfPnt2d( 1, aResult );
3534 theResultOnS2 = new TColgp_HArray1OfPnt2d( 1, aResult );
3535 theResultRadius = new TColStd_HArray1OfReal( 1, aResult );
3537 for ( Standard_Integer i = 1 ; i <= aResult; i++ ) {
3538 theResultOnS1->SetValue( i, aSeqResultS1.Value(i) );
3539 theResultOnS2->SetValue( i, aSeqResultS2.Value(i) );
3540 theResultRadius->SetValue( i, aSeqResultRad.Value(i) );
3546 // ------------------------------------------------------------------------------------------------
3547 // static function: AdjustByNeighbour
3549 // ------------------------------------------------------------------------------------------------
3550 gp_Pnt2d AdjustByNeighbour(const gp_Pnt2d& theaNeighbourPoint,
3551 const gp_Pnt2d& theOriginalPoint,
3552 Handle(GeomAdaptor_HSurface) theGASurface) {
3554 gp_Pnt2d ap1 = theaNeighbourPoint;
3555 gp_Pnt2d ap2 = theOriginalPoint;
3557 if ( theGASurface->IsUPeriodic() ) {
3558 Standard_Real aPeriod = theGASurface->UPeriod();
3559 gp_Pnt2d aPTest = ap2;
3560 Standard_Real aSqDistMin = 1.e+100;
3562 for ( Standard_Integer pIt = -1; pIt <= 1; pIt++) {
3563 aPTest.SetX( theOriginalPoint.X() + aPeriod * pIt );
3564 Standard_Real dd = ap1.SquareDistance( aPTest );
3566 if ( dd < aSqDistMin ) {
3572 if ( theGASurface->IsVPeriodic() ) {
3573 Standard_Real aPeriod = theGASurface->VPeriod();
3574 gp_Pnt2d aPTest = ap2;
3575 Standard_Real aSqDistMin = 1.e+100;
3577 for ( Standard_Integer pIt = -1; pIt <= 1; pIt++) {
3578 aPTest.SetY( theOriginalPoint.Y() + aPeriod * pIt );
3579 Standard_Real dd = ap1.SquareDistance( aPTest );
3581 if ( dd < aSqDistMin ) {
3590 // ------------------------------------------------------------------------------------------------
3591 //function: DecompositionOfWLine
3593 // ------------------------------------------------------------------------------------------------
3594 Standard_Boolean DecompositionOfWLine(const Handle(IntPatch_WLine)& theWLine,
3595 const Handle(GeomAdaptor_HSurface)& theSurface1,
3596 const Handle(GeomAdaptor_HSurface)& theSurface2,
3597 const TopoDS_Face& theFace1,
3598 const TopoDS_Face& theFace2,
3599 const GeomInt_LineConstructor& theLConstructor,
3600 const Standard_Boolean theAvoidLConstructor,
3601 IntPatch_SequenceOfLine& theNewLines,
3602 Standard_Real& theReachedTol3d,
3603 const Handle(IntTools_Context)& aContext)
3606 Standard_Boolean bRet, bAvoidLineConstructor;
3607 Standard_Integer aNbPnts, aNbParts;
3609 bRet=Standard_False;
3610 aNbPnts=theWLine->NbPnts();
3611 bAvoidLineConstructor=theAvoidLConstructor;
3616 if (!bAvoidLineConstructor) {
3617 aNbParts=theLConstructor.NbParts();
3623 Standard_Boolean bIsPrevPointOnBoundary, bIsPointOnBoundary, bIsCurrentPointOnBoundary;
3624 Standard_Integer nblines, pit, i, j;
3626 TColStd_Array1OfListOfInteger anArrayOfLines(1, aNbPnts);
3627 TColStd_Array1OfInteger anArrayOfLineType(1, aNbPnts);
3628 TColStd_ListOfInteger aListOfPointIndex;
3630 Handle(TColgp_HArray1OfPnt2d) aTanZoneS1;
3631 Handle(TColgp_HArray1OfPnt2d) aTanZoneS2;
3632 Handle(TColStd_HArray1OfReal) aTanZoneRadius;
3633 Standard_Integer aNbZone = ComputeTangentZones( theSurface1, theSurface2, theFace1, theFace2,
3634 aTanZoneS1, aTanZoneS2, aTanZoneRadius, aContext);
3638 aTol=Precision::Confusion();
3640 bIsPrevPointOnBoundary=Standard_False;
3641 bIsPointOnBoundary=Standard_False;
3646 for(pit = 1; pit <= aNbPnts; ++pit) {
3647 Standard_Boolean bIsOnFirstBoundary, isperiodic;
3648 Standard_Real aResolution, aPeriod, alowerboundary, aupperboundary, U, V;
3649 Standard_Real aParameter, anoffset, anAdjustPar;
3650 Standard_Real umin, umax, vmin, vmax;
3652 bIsCurrentPointOnBoundary = Standard_False;
3653 const IntSurf_PntOn2S& aPoint = theWLine->Point(pit);
3656 for(i = 0; i < 2; ++i) {
3657 Handle(GeomAdaptor_HSurface) aGASurface = (!i) ? theSurface1 : theSurface2;
3658 aGASurface->ChangeSurface().Surface()->Bounds(umin, umax, vmin, vmax);
3660 aPoint.ParametersOnS1(U, V);
3663 aPoint.ParametersOnS2(U, V);
3666 for(j = 0; j < 2; j++) {
3667 isperiodic = (!j) ? aGASurface->IsUPeriodic() : aGASurface->IsVPeriodic();
3673 aResolution=aGASurface->UResolution(aTol);
3674 aPeriod=aGASurface->UPeriod();
3675 alowerboundary=umin;
3676 aupperboundary=umax;
3680 aResolution=aGASurface->VResolution(aTol);
3681 aPeriod=aGASurface->VPeriod();
3682 alowerboundary=vmin;
3683 aupperboundary=vmax;
3687 GeomInt::AdjustPeriodic(aParameter,
3694 bIsOnFirstBoundary = Standard_True;// ?
3696 IsPointOnBoundary(anAdjustPar,
3700 bIsOnFirstBoundary);
3702 if(bIsPointOnBoundary) {
3703 bIsCurrentPointOnBoundary = Standard_True;
3707 // check if a point belong to a tangent zone. Begin
3708 Standard_Integer zIt = 0;
3709 for ( zIt = 1; zIt <= aNbZone; zIt++ ) {
3710 gp_Pnt2d aPZone = (i == 0) ? aTanZoneS1->Value(zIt) : aTanZoneS2->Value(zIt);
3711 Standard_Real aZoneRadius = aTanZoneRadius->Value(zIt);
3713 if ( IsInsideTanZone(gp_Pnt2d( U, V ), aPZone, aZoneRadius, aGASurface ) ) {
3714 // set boundary flag to split the curve by a tangent zone
3715 bIsPointOnBoundary = Standard_True;
3716 bIsCurrentPointOnBoundary = Standard_True;
3717 if ( theReachedTol3d < aZoneRadius ) {
3718 theReachedTol3d = aZoneRadius;
3724 }//for(j = 0; j < 2; j++) {
3726 if(bIsCurrentPointOnBoundary){
3729 }//for(i = 0; i < 2; ++i) {
3731 if((bIsCurrentPointOnBoundary != bIsPrevPointOnBoundary)) {
3732 if(!aListOfPointIndex.IsEmpty()) {
3734 anArrayOfLines.SetValue(nblines, aListOfPointIndex);
3735 anArrayOfLineType.SetValue(nblines, bIsPrevPointOnBoundary);
3736 aListOfPointIndex.Clear();
3738 bIsPrevPointOnBoundary = bIsCurrentPointOnBoundary;
3740 aListOfPointIndex.Append(pit);
3741 } //for(pit = 1; pit <= aNbPnts; ++pit) {
3743 if(!aListOfPointIndex.IsEmpty()) {
3745 anArrayOfLines.SetValue(nblines, aListOfPointIndex);
3746 anArrayOfLineType.SetValue(nblines, bIsPrevPointOnBoundary);
3747 aListOfPointIndex.Clear();
3751 return bRet; //Standard_False;
3755 // 2. Correct wlines.begin
3756 TColStd_Array1OfListOfInteger anArrayOfLineEnds(1, nblines);
3757 Handle(IntSurf_LineOn2S) aSeqOfPntOn2S = new IntSurf_LineOn2S();
3759 for(i = 1; i <= nblines; i++) {
3760 if(anArrayOfLineType.Value(i) != 0) {
3763 const TColStd_ListOfInteger& aListOfIndex = anArrayOfLines.Value(i);
3764 TColStd_ListOfInteger aListOfFLIndex;
3766 for(j = 0; j < 2; j++) {
3767 Standard_Integer aneighbourindex = (j == 0) ? (i - 1) : (i + 1);
3769 if((aneighbourindex < 1) || (aneighbourindex > nblines))
3772 if(anArrayOfLineType.Value(aneighbourindex) == 0)
3774 const TColStd_ListOfInteger& aNeighbour = anArrayOfLines.Value(aneighbourindex);
3775 Standard_Integer anIndex = (j == 0) ? aNeighbour.Last() : aNeighbour.First();
3776 const IntSurf_PntOn2S& aPoint = theWLine->Point(anIndex);
3778 IntSurf_PntOn2S aNewP = aPoint;
3779 if(aListOfIndex.Extent() < 2) {
3780 aSeqOfPntOn2S->Add(aNewP);
3781 aListOfFLIndex.Append(aSeqOfPntOn2S->NbPoints());
3785 Standard_Integer iFirst = aListOfIndex.First();
3786 Standard_Integer iLast = aListOfIndex.Last();
3788 for(Standard_Integer surfit = 0; surfit < 2; surfit++) {
3790 Handle(GeomAdaptor_HSurface) aGASurface = (surfit == 0) ? theSurface1 : theSurface2;
3791 Standard_Real umin=0., umax=0., vmin=0., vmax=0.;
3792 aGASurface->ChangeSurface().Surface()->Bounds(umin, umax, vmin, vmax);
3793 Standard_Real U=0., V=0.;
3796 aNewP.ParametersOnS1(U, V);
3798 aNewP.ParametersOnS2(U, V);
3799 Standard_Integer nbboundaries = 0;
3801 Standard_Boolean bIsNearBoundary = Standard_False;
3802 Standard_Integer aZoneIndex = 0;
3803 Standard_Integer bIsUBoundary = Standard_False; // use if nbboundaries == 1
3804 Standard_Integer bIsFirstBoundary = Standard_False; // use if nbboundaries == 1
3807 for(Standard_Integer parit = 0; parit < 2; parit++) {
3808 Standard_Boolean isperiodic = (parit == 0) ? aGASurface->IsUPeriodic() : aGASurface->IsVPeriodic();
3810 Standard_Real aResolution = (parit == 0) ? aGASurface->UResolution(aTol) : aGASurface->VResolution(aTol);
3811 Standard_Real alowerboundary = (parit == 0) ? umin : vmin;
3812 Standard_Real aupperboundary = (parit == 0) ? umax : vmax;
3814 Standard_Real aParameter = (parit == 0) ? U : V;
3815 Standard_Boolean bIsOnFirstBoundary = Standard_True;
3819 IsPointOnBoundary(aParameter, alowerboundary, aupperboundary, aResolution, bIsOnFirstBoundary);
3820 if(bIsPointOnBoundary) {
3821 bIsUBoundary = (parit == 0);
3822 bIsFirstBoundary = bIsOnFirstBoundary;
3827 Standard_Real aPeriod = (parit == 0) ? aGASurface->UPeriod() : aGASurface->VPeriod();
3828 Standard_Real anoffset, anAdjustPar;
3829 GeomInt::AdjustPeriodic(aParameter, alowerboundary, aupperboundary,
3830 aPeriod, anAdjustPar, anoffset);
3833 IsPointOnBoundary(anAdjustPar, alowerboundary, aupperboundary, aResolution, bIsOnFirstBoundary);
3834 if(bIsPointOnBoundary) {
3835 bIsUBoundary = (parit == 0);
3836 bIsFirstBoundary = bIsOnFirstBoundary;
3840 //check neighbourhood of boundary
3841 Standard_Real anEpsilon = aResolution * 100.;
3842 Standard_Real aPart = ( aupperboundary - alowerboundary ) * 0.1;
3843 anEpsilon = ( anEpsilon > aPart ) ? aPart : anEpsilon;
3845 bIsNearBoundary = IsPointOnBoundary(anAdjustPar, alowerboundary, aupperboundary,
3846 anEpsilon, bIsOnFirstBoundary);
3852 // check if a point belong to a tangent zone. Begin
3853 for ( Standard_Integer zIt = 1; zIt <= aNbZone; zIt++ ) {
3854 gp_Pnt2d aPZone = (surfit == 0) ? aTanZoneS1->Value(zIt) : aTanZoneS2->Value(zIt);
3855 Standard_Real aZoneRadius = aTanZoneRadius->Value(zIt);
3857 Standard_Integer aneighbourpointindex1 = (j == 0) ? iFirst : iLast;
3858 const IntSurf_PntOn2S& aNeighbourPoint = theWLine->Point(aneighbourpointindex1);
3859 Standard_Real nU1, nV1;
3862 aNeighbourPoint.ParametersOnS1(nU1, nV1);
3864 aNeighbourPoint.ParametersOnS2(nU1, nV1);
3865 gp_Pnt2d ap1(nU1, nV1);
3866 gp_Pnt2d ap2 = AdjustByNeighbour( ap1, gp_Pnt2d( U, V ), aGASurface );
3869 if ( IsInsideTanZone( ap2, aPZone, aZoneRadius, aGASurface ) ) {
3871 bIsNearBoundary = Standard_True;
3872 if ( theReachedTol3d < aZoneRadius ) {
3873 theReachedTol3d = aZoneRadius;
3877 // check if a point belong to a tangent zone. End
3878 Standard_Boolean bComputeLineEnd = Standard_False;
3880 if(nbboundaries == 2) {
3882 bComputeLineEnd = Standard_True;
3885 else if(nbboundaries == 1) {
3886 Standard_Boolean isperiodic = (bIsUBoundary) ? aGASurface->IsUPeriodic() : aGASurface->IsVPeriodic();
3889 Standard_Real alowerboundary = (bIsUBoundary) ? umin : vmin;
3890 Standard_Real aupperboundary = (bIsUBoundary) ? umax : vmax;
3891 Standard_Real aPeriod = (bIsUBoundary) ? aGASurface->UPeriod() : aGASurface->VPeriod();
3892 Standard_Real aParameter = (bIsUBoundary) ? U : V;
3893 Standard_Real anoffset, anAdjustPar;
3894 GeomInt::AdjustPeriodic(aParameter, alowerboundary, aupperboundary,
3895 aPeriod, anAdjustPar, anoffset);
3897 Standard_Real adist = (bIsFirstBoundary) ? fabs(anAdjustPar - alowerboundary) : fabs(anAdjustPar - aupperboundary);
3898 Standard_Real anotherPar = (bIsFirstBoundary) ? (aupperboundary - adist) : (alowerboundary + adist);
3899 anotherPar += anoffset;
3900 Standard_Integer aneighbourpointindex = (j == 0) ? aListOfIndex.First() : aListOfIndex.Last();
3901 const IntSurf_PntOn2S& aNeighbourPoint = theWLine->Point(aneighbourpointindex);
3902 Standard_Real nU1, nV1;
3905 aNeighbourPoint.ParametersOnS1(nU1, nV1);
3907 aNeighbourPoint.ParametersOnS2(nU1, nV1);
3909 Standard_Real adist1 = (bIsUBoundary) ? fabs(nU1 - U) : fabs(nV1 - V);
3910 Standard_Real adist2 = (bIsUBoundary) ? fabs(nU1 - anotherPar) : fabs(nV1 - anotherPar);
3911 bComputeLineEnd = Standard_True;
3912 Standard_Boolean bCheckAngle1 = Standard_False;
3913 Standard_Boolean bCheckAngle2 = Standard_False;
3915 Standard_Real anewU = (bIsUBoundary) ? anotherPar : U;
3916 Standard_Real anewV = (bIsUBoundary) ? V : anotherPar;
3918 if(((adist1 - adist2) > Precision::PConfusion()) &&
3919 (adist2 < (aPeriod / 4.))) {
3920 bCheckAngle1 = Standard_True;
3921 aNewVec = gp_Vec2d(gp_Pnt2d(nU1, nV1), gp_Pnt2d(anewU, anewV));
3923 if(aNewVec.SquareMagnitude() < gp::Resolution()) {
3924 aNewP.SetValue((surfit == 0), anewU, anewV);
3925 bCheckAngle1 = Standard_False;
3928 else if(adist1 < (aPeriod / 4.)) {
3929 bCheckAngle2 = Standard_True;
3930 aNewVec = gp_Vec2d(gp_Pnt2d(nU1, nV1), gp_Pnt2d(U, V));
3932 if(aNewVec.SquareMagnitude() < gp::Resolution()) {
3933 bCheckAngle2 = Standard_False;
3937 if(bCheckAngle1 || bCheckAngle2) {
3938 // assume there are at least two points in line (see "if" above)
3939 Standard_Integer anindexother = aneighbourpointindex;
3941 while((anindexother <= iLast) && (anindexother >= iFirst)) {
3942 anindexother = (j == 0) ? (anindexother + 1) : (anindexother - 1);
3943 const IntSurf_PntOn2S& aPrevNeighbourPoint = theWLine->Point(anindexother);
3944 Standard_Real nU2, nV2;
3947 aPrevNeighbourPoint.ParametersOnS1(nU2, nV2);
3949 aPrevNeighbourPoint.ParametersOnS2(nU2, nV2);
3950 gp_Vec2d aVecOld(gp_Pnt2d(nU2, nV2), gp_Pnt2d(nU1, nV1));
3952 if(aVecOld.SquareMagnitude() <= gp::Resolution()) {
3956 Standard_Real anAngle = aNewVec.Angle(aVecOld);
3958 if((fabs(anAngle) < (M_PI * 0.25)) && (aNewVec.Dot(aVecOld) > 0.)) {
3961 Standard_Real U1, U2, V1, V2;
3962 IntSurf_PntOn2S atmppoint = aNewP;
3963 atmppoint.SetValue((surfit == 0), anewU, anewV);
3964 atmppoint.Parameters(U1, V1, U2, V2);
3965 gp_Pnt P1 = theSurface1->Value(U1, V1);
3966 gp_Pnt P2 = theSurface2->Value(U2, V2);
3967 gp_Pnt P0 = aPoint.Value();
3969 if(P0.IsEqual(P1, aTol) &&
3970 P0.IsEqual(P2, aTol) &&
3971 P1.IsEqual(P2, aTol)) {
3972 bComputeLineEnd = Standard_False;
3973 aNewP.SetValue((surfit == 0), anewU, anewV);
3978 bComputeLineEnd = Standard_False;
3983 } // end while(anindexother...)
3987 else if ( bIsNearBoundary ) {
3988 bComputeLineEnd = Standard_True;
3991 if(bComputeLineEnd) {
3994 Standard_Boolean found = Standard_False;
3996 if ( bIsNearBoundary ) {
3997 // re-compute point near natural boundary or near tangent zone
3998 Standard_Real u1, v1, u2, v2;
3999 aNewP.Parameters( u1, v1, u2, v2 );
4001 anewpoint = gp_Pnt2d( u1, v1 );
4003 anewpoint = gp_Pnt2d( u2, v2 );
4005 Standard_Integer aneighbourpointindex1 = (j == 0) ? iFirst : iLast;
4006 const IntSurf_PntOn2S& aNeighbourPoint = theWLine->Point(aneighbourpointindex1);
4007 Standard_Real nU1, nV1;
4010 aNeighbourPoint.ParametersOnS1(nU1, nV1);
4012 aNeighbourPoint.ParametersOnS2(nU1, nV1);
4013 gp_Pnt2d ap1(nU1, nV1);
4018 // exclude point from a tangent zone
4019 anewpoint = AdjustByNeighbour( ap1, anewpoint, aGASurface );
4020 gp_Pnt2d aPZone = (surfit == 0) ? aTanZoneS1->Value(aZoneIndex) : aTanZoneS2->Value(aZoneIndex);
4021 Standard_Real aZoneRadius = aTanZoneRadius->Value(aZoneIndex);
4023 if ( FindPoint(ap1, anewpoint, umin, umax, vmin, vmax,
4024 aPZone, aZoneRadius, aGASurface, ap2) ) {
4026 found = Standard_True;
4029 else if ( aGASurface->IsUPeriodic() || aGASurface->IsVPeriodic() ) {
4030 // re-compute point near boundary if shifted on a period
4031 ap2 = AdjustByNeighbour( ap1, anewpoint, aGASurface );
4033 if ( ( ap2.X() < umin ) || ( ap2.X() > umax ) ||
4034 ( ap2.Y() < vmin ) || ( ap2.Y() > vmax ) ) {
4035 found = FindPoint(ap1, ap2, umin, umax, vmin, vmax, anewpoint);
4039 aNewP.SetValue( (surfit == 0), anewpoint.X(), anewpoint.Y() );
4045 Standard_Integer aneighbourpointindex1 = (j == 0) ? iFirst : iLast;
4046 const IntSurf_PntOn2S& aNeighbourPoint = theWLine->Point(aneighbourpointindex1);
4047 Standard_Real nU1, nV1;
4050 aNeighbourPoint.ParametersOnS1(nU1, nV1);
4052 aNeighbourPoint.ParametersOnS2(nU1, nV1);
4053 gp_Pnt2d ap1(nU1, nV1);
4054 gp_Pnt2d ap2(nU1, nV1);
4055 Standard_Integer aneighbourpointindex2 = aneighbourpointindex1;
4057 while((aneighbourpointindex2 <= iLast) && (aneighbourpointindex2 >= iFirst)) {
4058 aneighbourpointindex2 = (j == 0) ? (aneighbourpointindex2 + 1) : (aneighbourpointindex2 - 1);
4059 const IntSurf_PntOn2S& aPrevNeighbourPoint = theWLine->Point(aneighbourpointindex2);
4060 Standard_Real nU2, nV2;
4063 aPrevNeighbourPoint.ParametersOnS1(nU2, nV2);
4065 aPrevNeighbourPoint.ParametersOnS2(nU2, nV2);
4069 if(ap1.SquareDistance(ap2) > gp::Resolution()) {
4073 found = FindPoint(ap2, ap1, umin, umax, vmin, vmax, anewpoint);
4078 Standard_Real aCriteria = BRep_Tool::Tolerance(theFace1) + BRep_Tool::Tolerance(theFace2);
4079 GeomAPI_ProjectPointOnSurf& aProjector =
4080 (surfit == 0) ? aContext->ProjPS(theFace2) : aContext->ProjPS(theFace1);
4081 Handle(GeomAdaptor_HSurface) aSurface = (surfit == 0) ? theSurface1 : theSurface2;
4083 Handle(GeomAdaptor_HSurface) aSurfaceOther = (surfit == 0) ? theSurface2 : theSurface1;
4085 gp_Pnt aP3d = aSurface->Value(anewpoint.X(), anewpoint.Y());
4086 aProjector.Perform(aP3d);
4088 if(aProjector.IsDone()) {
4089 if(aProjector.LowerDistance() < aCriteria) {
4090 Standard_Real foundU = U, foundV = V;
4091 aProjector.LowerDistanceParameters(foundU, foundV);
4093 //Correction of projected coordinates. Begin
4094 //Note, it may be shifted on a period
4095 Standard_Integer aneindex1 = (j == 0) ? iFirst : iLast;
4096 const IntSurf_PntOn2S& aNeighbourPoint = theWLine->Point(aneindex1);
4097 Standard_Real nUn, nVn;
4100 aNeighbourPoint.ParametersOnS2(nUn, nVn);
4102 aNeighbourPoint.ParametersOnS1(nUn, nVn);
4103 gp_Pnt2d aNeighbour2d(nUn, nVn);
4104 gp_Pnt2d anAdjustedPoint = AdjustByNeighbour( aNeighbour2d, gp_Pnt2d(foundU, foundV), aSurfaceOther );
4105 foundU = anAdjustedPoint.X();
4106 foundV = anAdjustedPoint.Y();
4108 if ( ( anAdjustedPoint.X() < umin ) && ( anAdjustedPoint.X() > umax ) &&
4109 ( anAdjustedPoint.Y() < vmin ) && ( anAdjustedPoint.Y() > vmax ) ) {
4110 // attempt to roughly re-compute point
4111 foundU = ( foundU < umin ) ? umin : foundU;
4112 foundU = ( foundU > umax ) ? umax : foundU;
4113 foundV = ( foundV < vmin ) ? vmin : foundV;
4114 foundV = ( foundV > vmax ) ? vmax : foundV;
4116 GeomAPI_ProjectPointOnSurf& aProjector2 =
4117 (surfit == 0) ? aContext->ProjPS(theFace1) : aContext->ProjPS(theFace2);
4119 aP3d = aSurfaceOther->Value(foundU, foundV);
4120 aProjector2.Perform(aP3d);
4122 if(aProjector2.IsDone()) {
4123 if(aProjector2.LowerDistance() < aCriteria) {
4124 Standard_Real foundU2 = anewpoint.X(), foundV2 = anewpoint.Y();
4125 aProjector2.LowerDistanceParameters(foundU2, foundV2);
4126 anewpoint.SetX(foundU2);
4127 anewpoint.SetY(foundV2);
4131 //Correction of projected coordinates. End
4134 aNewP.SetValue(aP3d, anewpoint.X(), anewpoint.Y(), foundU, foundV);
4136 aNewP.SetValue(aP3d, foundU, foundV, anewpoint.X(), anewpoint.Y());
4142 aSeqOfPntOn2S->Add(aNewP);
4143 aListOfFLIndex.Append(aSeqOfPntOn2S->NbPoints());
4145 anArrayOfLineEnds.SetValue(i, aListOfFLIndex);
4147 // Correct wlines.end
4149 // Split wlines.begin
4150 Standard_Integer nbiter;
4153 if (!bAvoidLineConstructor) {
4154 nbiter=theLConstructor.NbParts();
4157 for(j = 1; j <= nbiter; ++j) {
4158 Standard_Real fprm, lprm;
4159 Standard_Integer ifprm, ilprm;
4161 if(bAvoidLineConstructor) {
4163 ilprm = theWLine->NbPnts();
4166 theLConstructor.Part(j, fprm, lprm);
4167 ifprm = (Standard_Integer)fprm;
4168 ilprm = (Standard_Integer)lprm;
4171 Handle(IntSurf_LineOn2S) aLineOn2S = new IntSurf_LineOn2S();
4173 for(i = 1; i <= nblines; i++) {
4174 if(anArrayOfLineType.Value(i) != 0) {
4177 const TColStd_ListOfInteger& aListOfIndex = anArrayOfLines.Value(i);
4178 const TColStd_ListOfInteger& aListOfFLIndex = anArrayOfLineEnds.Value(i);
4179 Standard_Boolean bhasfirstpoint = (aListOfFLIndex.Extent() == 2);
4180 Standard_Boolean bhaslastpoint = (aListOfFLIndex.Extent() == 2);
4182 if(!bhasfirstpoint && !aListOfFLIndex.IsEmpty()) {
4183 bhasfirstpoint = (i != 1);
4186 if(!bhaslastpoint && !aListOfFLIndex.IsEmpty()) {
4187 bhaslastpoint = (i != nblines);
4190 Standard_Integer iFirst = aListOfIndex.First();
4191 Standard_Integer iLast = aListOfIndex.Last();
4192 Standard_Boolean bIsFirstInside = ((ifprm >= iFirst) && (ifprm <= iLast));
4193 Standard_Boolean bIsLastInside = ((ilprm >= iFirst) && (ilprm <= iLast));
4195 if(!bIsFirstInside && !bIsLastInside) {
4196 if((ifprm < iFirst) && (ilprm > iLast)) {
4197 // append whole line, and boundaries if neccesary
4198 if(bhasfirstpoint) {
4199 pit = aListOfFLIndex.First();
4200 const IntSurf_PntOn2S& aP = aSeqOfPntOn2S->Value(pit);
4203 TColStd_ListIteratorOfListOfInteger anIt(aListOfIndex);
4205 for(; anIt.More(); anIt.Next()) {
4207 const IntSurf_PntOn2S& aP = theWLine->Point(pit);
4212 pit = aListOfFLIndex.Last();
4213 const IntSurf_PntOn2S& aP = aSeqOfPntOn2S->Value(pit);
4217 // check end of split line (end is almost always)
4218 Standard_Integer aneighbour = i + 1;
4219 Standard_Boolean bIsEndOfLine = Standard_True;
4221 if(aneighbour <= nblines) {
4222 const TColStd_ListOfInteger& aListOfNeighbourIndex = anArrayOfLines.Value(aneighbour);
4224 if((anArrayOfLineType.Value(aneighbour) != 0) &&
4225 (aListOfNeighbourIndex.IsEmpty())) {
4226 bIsEndOfLine = Standard_False;
4231 if(aLineOn2S->NbPoints() > 1) {
4232 Handle(IntPatch_WLine) aNewWLine =
4233 new IntPatch_WLine(aLineOn2S, Standard_False);
4234 theNewLines.Append(aNewWLine);
4236 aLineOn2S = new IntSurf_LineOn2S();
4241 // end if(!bIsFirstInside && !bIsLastInside)
4243 if(bIsFirstInside && bIsLastInside) {
4244 // append inside points between ifprm and ilprm
4245 TColStd_ListIteratorOfListOfInteger anIt(aListOfIndex);
4247 for(; anIt.More(); anIt.Next()) {
4249 if((pit < ifprm) || (pit > ilprm))
4251 const IntSurf_PntOn2S& aP = theWLine->Point(pit);
4257 if(bIsFirstInside) {
4258 // append points from ifprm to last point + boundary point
4259 TColStd_ListIteratorOfListOfInteger anIt(aListOfIndex);
4261 for(; anIt.More(); anIt.Next()) {
4265 const IntSurf_PntOn2S& aP = theWLine->Point(pit);
4270 pit = aListOfFLIndex.Last();
4271 const IntSurf_PntOn2S& aP = aSeqOfPntOn2S->Value(pit);
4274 // check end of split line (end is almost always)
4275 Standard_Integer aneighbour = i + 1;
4276 Standard_Boolean bIsEndOfLine = Standard_True;
4278 if(aneighbour <= nblines) {
4279 const TColStd_ListOfInteger& aListOfNeighbourIndex = anArrayOfLines.Value(aneighbour);
4281 if((anArrayOfLineType.Value(aneighbour) != 0) &&
4282 (aListOfNeighbourIndex.IsEmpty())) {
4283 bIsEndOfLine = Standard_False;
4288 if(aLineOn2S->NbPoints() > 1) {
4289 Handle(IntPatch_WLine) aNewWLine =
4290 new IntPatch_WLine(aLineOn2S, Standard_False);
4291 theNewLines.Append(aNewWLine);
4293 aLineOn2S = new IntSurf_LineOn2S();
4296 // end if(bIsFirstInside)
4299 // append points from first boundary point to ilprm
4300 if(bhasfirstpoint) {
4301 pit = aListOfFLIndex.First();
4302 const IntSurf_PntOn2S& aP = aSeqOfPntOn2S->Value(pit);
4305 TColStd_ListIteratorOfListOfInteger anIt(aListOfIndex);
4307 for(; anIt.More(); anIt.Next()) {
4311 const IntSurf_PntOn2S& aP = theWLine->Point(pit);
4315 //end if(bIsLastInside)
4319 if(aLineOn2S->NbPoints() > 1) {
4320 Handle(IntPatch_WLine) aNewWLine =
4321 new IntPatch_WLine(aLineOn2S, Standard_False);
4322 theNewLines.Append(aNewWLine);
4327 return Standard_True;
4330 // ------------------------------------------------------------------------------------------------
4331 // static function: ParameterOutOfBoundary
4332 // purpose: Computes a new parameter for given curve. The corresponding 2d points
4333 // does not lay on any boundary of given faces
4334 // ------------------------------------------------------------------------------------------------
4335 Standard_Boolean ParameterOutOfBoundary(const Standard_Real theParameter,
4336 const Handle(Geom_Curve)& theCurve,
4337 const TopoDS_Face& theFace1,
4338 const TopoDS_Face& theFace2,
4339 const Standard_Real theOtherParameter,
4340 const Standard_Boolean bIncreasePar,
4341 Standard_Real& theNewParameter,
4342 const Handle(IntTools_Context)& aContext)
4344 Standard_Boolean bIsComputed = Standard_False;
4345 theNewParameter = theParameter;
4347 Standard_Real acurpar = theParameter;
4348 TopAbs_State aState = TopAbs_ON;
4349 Standard_Integer iter = 0;
4350 Standard_Real asumtol = BRep_Tool::Tolerance(theFace1) + BRep_Tool::Tolerance(theFace2);
4351 Standard_Real adelta = asumtol * 0.1;
4352 adelta = (adelta < Precision::Confusion()) ? Precision::Confusion() : adelta;
4353 Handle(Geom_Surface) aSurf1 = BRep_Tool::Surface(theFace1);
4354 Handle(Geom_Surface) aSurf2 = BRep_Tool::Surface(theFace2);
4356 Standard_Real u1, u2, v1, v2;
4358 GeomAPI_ProjectPointOnSurf aPrj1;
4359 aSurf1->Bounds(u1, u2, v1, v2);
4360 aPrj1.Init(aSurf1, u1, u2, v1, v2);
4362 GeomAPI_ProjectPointOnSurf aPrj2;
4363 aSurf2->Bounds(u1, u2, v1, v2);
4364 aPrj2.Init(aSurf2, u1, u2, v1, v2);
4366 while(aState == TopAbs_ON) {
4371 gp_Pnt aPCurrent = theCurve->Value(acurpar);
4372 aPrj1.Perform(aPCurrent);
4373 Standard_Real U=0., V=0.;
4375 if(aPrj1.IsDone()) {
4376 aPrj1.LowerDistanceParameters(U, V);
4377 aState = aContext->StatePointFace(theFace1, gp_Pnt2d(U, V));
4380 if(aState != TopAbs_ON) {
4381 aPrj2.Perform(aPCurrent);
4383 if(aPrj2.IsDone()) {
4384 aPrj2.LowerDistanceParameters(U, V);
4385 aState = aContext->StatePointFace(theFace2, gp_Pnt2d(U, V));
4396 theNewParameter = acurpar;
4397 bIsComputed = Standard_True;
4400 if(acurpar >= theOtherParameter)
4401 theNewParameter = theOtherParameter;
4404 if(acurpar <= theOtherParameter)
4405 theNewParameter = theOtherParameter;
4411 //=======================================================================
4412 //function : IsCurveValid
4414 //=======================================================================
4415 Standard_Boolean IsCurveValid(Handle(Geom2d_Curve)& thePCurve)
4417 if(thePCurve.IsNull())
4418 return Standard_False;
4420 Standard_Real tolint = 1.e-10;
4421 Geom2dAdaptor_Curve PCA;
4422 IntRes2d_Domain PCD;
4423 Geom2dInt_GInter PCI;
4425 Standard_Real pf = 0., pl = 0.;
4426 gp_Pnt2d pntf, pntl;
4428 if(!thePCurve->IsClosed() && !thePCurve->IsPeriodic()) {
4429 pf = thePCurve->FirstParameter();
4430 pl = thePCurve->LastParameter();
4431 pntf = thePCurve->Value(pf);
4432 pntl = thePCurve->Value(pl);
4433 PCA.Load(thePCurve);
4434 if(!PCA.IsPeriodic()) {
4435 if(PCA.FirstParameter() > pf) pf = PCA.FirstParameter();
4436 if(PCA.LastParameter() < pl) pl = PCA.LastParameter();
4438 PCD.SetValues(pntf,pf,tolint,pntl,pl,tolint);
4439 PCI.Perform(PCA,PCD,tolint,tolint);
4441 if(PCI.NbPoints() > 0) {
4442 return Standard_False;
4446 return Standard_True;
4449 //=======================================================================
4450 //static function : ApproxWithPCurves
4451 //purpose : for bug 20964 only
4452 //=======================================================================
4453 Standard_Boolean ApproxWithPCurves(const gp_Cylinder& theCyl,
4454 const gp_Sphere& theSph)
4456 Standard_Boolean bRes = Standard_True;
4457 Standard_Real R1 = theCyl.Radius(), R2 = theSph.Radius();
4460 Standard_Real aD2, aRc2, aEps;
4466 const gp_Ax3& aAx3Sph=theSph.Position();
4467 const gp_Pnt& aLocSph=aAx3Sph.Location();
4468 const gp_Dir& aDirSph=aAx3Sph.Direction();
4470 const gp_Ax1& aAx1Cyl=theCyl.Axis();
4471 gp_Lin aLinCyl(aAx1Cyl);
4473 aApexSph.SetXYZ(aLocSph.XYZ()+R2*aDirSph.XYZ());
4474 aD2=aLinCyl.SquareDistance(aApexSph);
4475 if (fabs(aD2-aRc2)<aEps) {
4479 aApexSph.SetXYZ(aLocSph.XYZ()-R2*aDirSph.XYZ());
4480 aD2=aLinCyl.SquareDistance(aApexSph);
4481 if (fabs(aD2-aRc2)<aEps) {
4490 gp_Lin anCylAx(theCyl.Axis());
4492 Standard_Real aDist = anCylAx.Distance(theSph.Location());
4493 Standard_Real aDRel = Abs(aDist - R1)/R2;
4495 if(aDRel > .2) return bRes;
4497 Standard_Real par = ElCLib::Parameter(anCylAx, theSph.Location());
4498 gp_Pnt aP = ElCLib::Value(par, anCylAx);
4499 gp_Vec aV(aP, theSph.Location());
4501 Standard_Real dd = aV.Dot(theSph.Position().XDirection());
4503 if(aDist < R1 && dd > 0.) return Standard_False;
4504 if(aDist > R1 && dd < 0.) return Standard_False;
4509 //=======================================================================
4510 //function : PerformPlanes
4512 //=======================================================================
4513 void PerformPlanes(const Handle(GeomAdaptor_HSurface)& theS1,
4514 const Handle(GeomAdaptor_HSurface)& theS2,
4515 const Standard_Real TolAng,
4516 const Standard_Real TolTang,
4517 const Standard_Boolean theApprox1,
4518 const Standard_Boolean theApprox2,
4519 IntTools_SequenceOfCurves& theSeqOfCurve,
4520 Standard_Boolean& theTangentFaces)
4523 gp_Pln aPln1 = theS1->Surface().Plane();
4524 gp_Pln aPln2 = theS2->Surface().Plane();
4526 IntAna_QuadQuadGeo aPlnInter(aPln1, aPln2, TolAng, TolTang);
4528 if(!aPlnInter.IsDone()) {
4529 theTangentFaces = Standard_False;
4533 IntAna_ResultType aResType = aPlnInter.TypeInter();
4535 if(aResType == IntAna_Same) {
4536 theTangentFaces = Standard_True;
4540 theTangentFaces = Standard_False;
4542 if(aResType == IntAna_Empty) {
4546 gp_Lin aLin = aPlnInter.Line(1);
4548 ProjLib_Plane aProj;
4551 aProj.Project(aLin);
4552 gp_Lin2d aLin2d1 = aProj.Line();
4555 aProj.Project(aLin);
4556 gp_Lin2d aLin2d2 = aProj.Line();
4558 //classify line2d1 relatively first plane
4559 Standard_Real P11, P12;
4560 Standard_Boolean IsCrossed = ClassifyLin2d(theS1, aLin2d1, TolTang, P11, P12);
4561 if(!IsCrossed) return;
4562 //classify line2d2 relatively second plane
4563 Standard_Real P21, P22;
4564 IsCrossed = ClassifyLin2d(theS2, aLin2d2, TolTang, P21, P22);
4565 if(!IsCrossed) return;
4567 //Analysis of parametric intervals: must have common part
4569 if(P21 >= P12) return;
4570 if(P22 <= P11) return;
4572 Standard_Real pmin, pmax;
4573 pmin = Max(P11, P21);
4574 pmax = Min(P12, P22);
4576 if(pmax - pmin <= TolTang) return;
4578 Handle(Geom_Line) aGLin = new Geom_Line(aLin);
4580 IntTools_Curve aCurve;
4581 Handle(Geom_TrimmedCurve) aGTLin = new Geom_TrimmedCurve(aGLin, pmin, pmax);
4583 aCurve.SetCurve(aGTLin);
4586 Handle(Geom2d_Line) C2d = new Geom2d_Line(aLin2d1);
4587 aCurve.SetFirstCurve2d(new Geom2d_TrimmedCurve(C2d, pmin, pmax));
4590 Handle(Geom2d_Curve) H1;
4591 aCurve.SetFirstCurve2d(H1);
4594 Handle(Geom2d_Line) C2d = new Geom2d_Line(aLin2d2);
4595 aCurve.SetSecondCurve2d(new Geom2d_TrimmedCurve(C2d, pmin, pmax));
4598 Handle(Geom2d_Curve) H1;
4599 aCurve.SetFirstCurve2d(H1);
4602 theSeqOfCurve.Append(aCurve);
4606 //=======================================================================
4607 //function : ClassifyLin2d
4609 //=======================================================================
4610 static inline Standard_Boolean INTER(const Standard_Real d1,
4611 const Standard_Real d2,
4612 const Standard_Real tol)
4614 return (d1 > tol && d2 < -tol) ||
4615 (d1 < -tol && d2 > tol) ||
4616 ((d1 <= tol && d1 >= -tol) && (d2 > tol || d2 < -tol)) ||
4617 ((d2 <= tol && d2 >= -tol) && (d1 > tol || d1 < -tol));
4619 static inline Standard_Boolean COINC(const Standard_Real d1,
4620 const Standard_Real d2,
4621 const Standard_Real tol)
4623 return (d1 <= tol && d1 >= -tol) && (d2 <= tol && d2 >= -tol);
4625 Standard_Boolean ClassifyLin2d(const Handle(GeomAdaptor_HSurface)& theS,
4626 const gp_Lin2d& theLin2d,
4627 const Standard_Real theTol,
4628 Standard_Real& theP1,
4629 Standard_Real& theP2)
4632 Standard_Real xmin, xmax, ymin, ymax, d1, d2, A, B, C;
4633 Standard_Real par[2];
4634 Standard_Integer nbi = 0;
4636 xmin = theS->Surface().FirstUParameter();
4637 xmax = theS->Surface().LastUParameter();
4638 ymin = theS->Surface().FirstVParameter();
4639 ymax = theS->Surface().LastVParameter();
4641 theLin2d.Coefficients(A, B, C);
4643 //xmin, ymin <-> xmin, ymax
4644 d1 = A*xmin + B*ymin + C;
4645 d2 = A*xmin + B*ymax + C;
4647 if(INTER(d1, d2, theTol)) {
4648 //Intersection with boundary
4649 Standard_Real y = -(C + A*xmin)/B;
4650 par[nbi] = ElCLib::Parameter(theLin2d, gp_Pnt2d(xmin, y));
4653 else if (COINC(d1, d2, theTol)) {
4654 //Coincidence with boundary
4655 par[0] = ElCLib::Parameter(theLin2d, gp_Pnt2d(xmin, ymin));
4656 par[1] = ElCLib::Parameter(theLin2d, gp_Pnt2d(xmin, ymax));
4662 if(fabs(par[0]-par[1]) > theTol) {
4663 theP1 = Min(par[0], par[1]);
4664 theP2 = Max(par[0], par[1]);
4665 return Standard_True;
4667 else return Standard_False;
4671 //xmin, ymax <-> xmax, ymax
4673 d2 = A*xmax + B*ymax + C;
4675 if(d1 > theTol || d1 < -theTol) {//to avoid checking of
4676 //coincidence with the same point
4677 if(INTER(d1, d2, theTol)) {
4678 Standard_Real x = -(C + B*ymax)/A;
4679 par[nbi] = ElCLib::Parameter(theLin2d, gp_Pnt2d(x, ymax));
4682 else if (COINC(d1, d2, theTol)) {
4683 par[0] = ElCLib::Parameter(theLin2d, gp_Pnt2d(xmin, ymax));
4684 par[1] = ElCLib::Parameter(theLin2d, gp_Pnt2d(xmax, ymax));
4691 if(fabs(par[0]-par[1]) > theTol) {
4692 theP1 = Min(par[0], par[1]);
4693 theP2 = Max(par[0], par[1]);
4694 return Standard_True;
4696 else return Standard_False;
4700 //xmax, ymax <-> xmax, ymin
4702 d2 = A*xmax + B*ymin + C;
4704 if(d1 > theTol || d1 < -theTol) {
4705 if(INTER(d1, d2, theTol)) {
4706 Standard_Real y = -(C + A*xmax)/B;
4707 par[nbi] = ElCLib::Parameter(theLin2d, gp_Pnt2d(xmax, y));
4710 else if (COINC(d1, d2, theTol)) {
4711 par[0] = ElCLib::Parameter(theLin2d, gp_Pnt2d(xmax, ymax));
4712 par[1] = ElCLib::Parameter(theLin2d, gp_Pnt2d(xmax, ymin));
4718 if(fabs(par[0]-par[1]) > theTol) {
4719 theP1 = Min(par[0], par[1]);
4720 theP2 = Max(par[0], par[1]);
4721 return Standard_True;
4723 else return Standard_False;
4726 //xmax, ymin <-> xmin, ymin
4728 d2 = A*xmin + B*ymin + C;
4730 if(d1 > theTol || d1 < -theTol) {
4731 if(INTER(d1, d2, theTol)) {
4732 Standard_Real x = -(C + B*ymin)/A;
4733 par[nbi] = ElCLib::Parameter(theLin2d, gp_Pnt2d(x, ymin));
4736 else if (COINC(d1, d2, theTol)) {
4737 par[0] = ElCLib::Parameter(theLin2d, gp_Pnt2d(xmax, ymin));
4738 par[1] = ElCLib::Parameter(theLin2d, gp_Pnt2d(xmin, ymin));
4744 if(fabs(par[0]-par[1]) > theTol) {
4745 theP1 = Min(par[0], par[1]);
4746 theP2 = Max(par[0], par[1]);
4747 return Standard_True;
4749 else return Standard_False;
4752 return Standard_False;
4756 //=======================================================================
4757 //function : ApproxParameters
4759 //=======================================================================
4760 void ApproxParameters(const Handle(GeomAdaptor_HSurface)& aHS1,
4761 const Handle(GeomAdaptor_HSurface)& aHS2,
4762 Standard_Integer& iDegMin,
4763 Standard_Integer& iDegMax,
4764 Standard_Integer& iNbIter)
4767 GeomAbs_SurfaceType aTS1, aTS2;
4774 aTS1=aHS1->Surface().GetType();
4775 aTS2=aHS2->Surface().GetType();
4778 if ((aTS1==GeomAbs_Cylinder && aTS2==GeomAbs_Torus) ||
4779 (aTS2==GeomAbs_Cylinder && aTS1==GeomAbs_Torus)) {
4780 Standard_Real aRC, aRT, dR, aPC;
4781 gp_Cylinder aCylinder;
4784 aPC=Precision::Confusion();
4786 aCylinder=(aTS1==GeomAbs_Cylinder)? aHS1->Surface().Cylinder() : aHS2->Surface().Cylinder();
4787 aTorus=(aTS1==GeomAbs_Torus)? aHS1->Surface().Torus() : aHS2->Surface().Torus();
4789 aRC=aCylinder.Radius();
4790 aRT=aTorus.MinorRadius();
4800 if (aTS1==GeomAbs_Cylinder && aTS2==GeomAbs_Cylinder) {
4804 //=======================================================================
4805 //function : Tolerances
4807 //=======================================================================
4808 void Tolerances(const Handle(GeomAdaptor_HSurface)& aHS1,
4809 const Handle(GeomAdaptor_HSurface)& aHS2,
4810 Standard_Real& aTolTang)
4812 GeomAbs_SurfaceType aTS1, aTS2;
4814 aTS1=aHS1->Surface().GetType();
4815 aTS2=aHS2->Surface().GetType();
4818 if ((aTS1==GeomAbs_Cylinder && aTS2==GeomAbs_Torus) ||
4819 (aTS2==GeomAbs_Cylinder && aTS1==GeomAbs_Torus)) {
4820 Standard_Real aRC, aRT, dR, aPC;
4821 gp_Cylinder aCylinder;
4824 aPC=Precision::Confusion();
4826 aCylinder=(aTS1==GeomAbs_Cylinder)? aHS1->Surface().Cylinder() : aHS2->Surface().Cylinder();
4827 aTorus=(aTS1==GeomAbs_Torus)? aHS1->Surface().Torus() : aHS2->Surface().Torus();
4829 aRC=aCylinder.Radius();
4830 aRT=aTorus.MinorRadius();
4837 aTolTang=0.1*aTolTang;
4841 //=======================================================================
4842 //function : SortTypes
4844 //=======================================================================
4845 Standard_Boolean SortTypes(const GeomAbs_SurfaceType aType1,
4846 const GeomAbs_SurfaceType aType2)
4848 Standard_Boolean bRet;
4849 Standard_Integer aI1, aI2;
4851 bRet=Standard_False;
4853 aI1=IndexType(aType1);
4854 aI2=IndexType(aType2);
4860 //=======================================================================
4861 //function : IndexType
4863 //=======================================================================
4864 Standard_Integer IndexType(const GeomAbs_SurfaceType aType)
4866 Standard_Integer aIndex;
4870 if (aType==GeomAbs_Plane) {
4873 else if (aType==GeomAbs_Cylinder) {
4876 else if (aType==GeomAbs_Cone) {
4879 else if (aType==GeomAbs_Sphere) {
4882 else if (aType==GeomAbs_Torus) {
4885 else if (aType==GeomAbs_BezierSurface) {
4888 else if (aType==GeomAbs_BSplineSurface) {
4891 else if (aType==GeomAbs_SurfaceOfRevolution) {
4894 else if (aType==GeomAbs_SurfaceOfExtrusion) {
4897 else if (aType==GeomAbs_OffsetSurface) {
4900 else if (aType==GeomAbs_OtherSurface) {
4905 #ifdef OCCT_DEBUG_DUMPWLINE
4906 //=======================================================================
4907 //function : DumpWLine
4909 //=======================================================================
4910 void DumpWLine(const Handle(IntPatch_WLine)& aWLine)
4912 Standard_Integer i, aNbPnts;
4913 Standard_Real aX, aY, aZ, aU1, aV1, aU2, aV2;
4915 printf(" *WLine\n");
4916 aNbPnts=aWLine->NbPnts();
4917 for (i=1; i<=aNbPnts; ++i) {
4918 const IntSurf_PntOn2S aPntOn2S=aWLine->Point(i);
4919 const gp_Pnt& aP3D=aPntOn2S.Value();
4920 aP3D.Coord(aX, aY, aZ);
4921 aPntOn2S.Parameters(aU1, aV1, aU2, aV2);
4923 printf("point p_%d %lf %lf %lf\n", i, aX, aY, aZ);
4924 //printf("point p_%d %20.15lf %20.15lf %20.15lf %20.15lf %20.15lf %20.15lf %20.15lf\n",
4925 // i, aX, aY, aZ, aU1, aV1, aU2, aV2);
4929 //=======================================================================
4930 //function : RefineVector
4932 //=======================================================================
4933 void RefineVector(gp_Vec2d& aV2D)
4935 Standard_Integer k,m;
4936 Standard_Real aC[2], aEps, aR1, aR2, aNum;
4942 aV2D.Coord(aC[0], aC[1]);
4944 for (k=0; k<2; ++k) {
4947 if (aNum>aR1 && aNum<aR2) {
4958 aV2D.SetCoord(aC[0], aC[1]);
4961 //=======================================================================
4962 // Function : FindMaxDistance
4964 //=======================================================================
4965 Standard_Real FindMaxDistance(const Handle(Geom_Curve)& theCurve,
4966 const Standard_Real theFirst,
4967 const Standard_Real theLast,
4968 const TopoDS_Face& theFace,
4969 const Handle(IntTools_Context)& theContext)
4971 Standard_Integer aNbS;
4972 Standard_Real aT1, aT2, aDt, aD, aDMax, anEps;
4975 aDt = (theLast - theFirst) / aNbS;
4977 anEps = 1.e-4 * aDt;
4979 GeomAPI_ProjectPointOnSurf& aProjPS = theContext->ProjPS(theFace);
4985 if (aT2 > theLast) {
4989 aD = FindMaxDistance(theCurve, aT1, aT2, aProjPS, anEps);
4998 //=======================================================================
4999 // Function : FindMaxDistance
5001 //=======================================================================
5002 Standard_Real FindMaxDistance(const Handle(Geom_Curve)& theC,
5003 const Standard_Real theFirst,
5004 const Standard_Real theLast,
5005 GeomAPI_ProjectPointOnSurf& theProjPS,
5006 const Standard_Real theEps)
5008 Standard_Real aA, aB, aCf, aX, aX1, aX2, aF1, aF2, aF;
5010 aCf = 0.61803398874989484820458683436564;//(sqrt(5.)-1)/2.;
5014 aX1 = aB - aCf * (aB - aA);
5015 aF1 = MaxDistance(theC, aX1, theProjPS);
5016 aX2 = aA + aCf * (aB - aA);
5017 aF2 = MaxDistance(theC, aX2, theProjPS);
5020 if ((aB - aA) < theEps) {
5028 aX1 = aB - aCf * (aB - aA);
5029 aF1 = MaxDistance(theC, aX1, theProjPS);
5035 aX2 = aA + aCf * (aB - aA);
5036 aF2 = MaxDistance(theC, aX2, theProjPS);
5040 aX = 0.5 * (aA + aB);
5041 aF = MaxDistance(theC, aX, theProjPS);
5054 //=======================================================================
5055 // Function : MaxDistance
5057 //=======================================================================
5058 Standard_Real MaxDistance(const Handle(Geom_Curve)& theC,
5059 const Standard_Real aT,
5060 GeomAPI_ProjectPointOnSurf& theProjPS)
5066 theProjPS.Perform(aP);
5067 aD = theProjPS.NbPoints() ? theProjPS.LowerDistance() : 0.;
5072 //=======================================================================
5073 //function : CheckPCurve
5074 //purpose : Checks if points of the pcurve are out of the face bounds.
5075 //=======================================================================
5076 Standard_Boolean CheckPCurve(const Handle(Geom2d_Curve)& aPC,
5077 const TopoDS_Face& aFace)
5079 const Standard_Integer NPoints = 23;
5081 Standard_Real umin,umax,vmin,vmax;
5083 BRepTools::UVBounds(aFace, umin, umax, vmin, vmax);
5084 Standard_Real tolU = Max ((umax-umin)*0.01, Precision::Confusion());
5085 Standard_Real tolV = Max ((vmax-vmin)*0.01, Precision::Confusion());
5086 Standard_Real fp = aPC->FirstParameter();
5087 Standard_Real lp = aPC->LastParameter();
5090 // adjust domain for periodic surfaces
5091 TopLoc_Location aLoc;
5092 Handle(Geom_Surface) aSurf = BRep_Tool::Surface(aFace, aLoc);
5093 if (aSurf->IsKind(STANDARD_TYPE(Geom_RectangularTrimmedSurface))) {
5094 aSurf = (Handle(Geom_RectangularTrimmedSurface)::DownCast(aSurf))->BasisSurface();
5096 gp_Pnt2d pnt = aPC->Value((fp+lp)/2);
5100 if (aSurf->IsUPeriodic()) {
5101 Standard_Real aPer = aSurf->UPeriod();
5102 Standard_Integer nshift = (Standard_Integer) ((u-umin)/aPer);
5103 if (u < umin+aPer*nshift) nshift--;
5104 umin += aPer*nshift;
5105 umax += aPer*nshift;
5107 if (aSurf->IsVPeriodic()) {
5108 Standard_Real aPer = aSurf->VPeriod();
5109 Standard_Integer nshift = (Standard_Integer) ((v-vmin)/aPer);
5110 if (v < vmin+aPer*nshift) nshift--;
5111 vmin += aPer*nshift;
5112 vmax += aPer*nshift;
5115 //--------------------------------------------------------
5116 Standard_Boolean bRet;
5117 Standard_Integer j, aNbIntervals;
5118 Standard_Real aT, dT;
5121 Geom2dAdaptor_Curve aGAC(aPC);
5122 aNbIntervals=aGAC.NbIntervals(GeomAbs_CN);
5124 TColStd_Array1OfReal aTI(1, aNbIntervals+1);
5125 aGAC.Intervals(aTI,GeomAbs_CN);
5127 bRet=Standard_False;
5129 aT=aGAC.FirstParameter();
5130 for (j=1; j<=aNbIntervals; ++j) {
5131 dT=(aTI(j+1)-aTI(j))/NPoints;
5133 for (i=1; i<NPoints; i++) {
5137 if (umin-u > tolU || u-umax > tolU ||
5138 vmin-v > tolV || v-vmax > tolV) {
5145 //=======================================================================
5146 //function : CorrectPlaneBoundaries
5148 //=======================================================================
5149 void CorrectPlaneBoundaries(Standard_Real& aUmin,
5150 Standard_Real& aUmax,
5151 Standard_Real& aVmin,
5152 Standard_Real& aVmax)
5154 if (!(Precision::IsInfinite(aUmin) ||
5155 Precision::IsInfinite(aUmax))) {
5158 dU=0.1*(aUmax-aUmin);
5162 if (!(Precision::IsInfinite(aVmin) ||
5163 Precision::IsInfinite(aVmax))) {
5166 dV=0.1*(aVmax-aVmin);