0023405: BOP common produces one face instead of a solid
[occt.git] / src / IntTools / IntTools_FaceFace.cxx
... / ...
CommitLineData
1// Created on: 2000-11-23
2// Created by: Michael KLOKOV
3// Copyright (c) 2000-2012 OPEN CASCADE SAS
4//
5// The content of this file is subject to the Open CASCADE Technology Public
6// License Version 6.5 (the "License"). You may not use the content of this file
7// except in compliance with the License. Please obtain a copy of the License
8// at http://www.opencascade.org and read it completely before using this file.
9//
10// The Initial Developer of the Original Code is Open CASCADE S.A.S., having its
11// main offices at: 1, place des Freres Montgolfier, 78280 Guyancourt, France.
12//
13// The Original Code and all software distributed under the License is
14// distributed on an "AS IS" basis, without warranty of any kind, and the
15// Initial Developer hereby disclaims all such warranties, including without
16// limitation, any warranties of merchantability, fitness for a particular
17// purpose or non-infringement. Please see the License for the specific terms
18// and conditions governing the rights and limitations under the License.
19
20
21
22#include <IntTools_FaceFace.ixx>
23
24#include <Precision.hxx>
25
26#include <TColStd_HArray1OfReal.hxx>
27#include <TColStd_Array1OfReal.hxx>
28#include <TColStd_Array1OfInteger.hxx>
29#include <TColStd_SequenceOfReal.hxx>
30#include <TColStd_ListOfInteger.hxx>
31#include <TColStd_ListIteratorOfListOfInteger.hxx>
32#include <TColStd_Array1OfListOfInteger.hxx>
33
34#include <gp_Lin2d.hxx>
35#include <gp_Ax22d.hxx>
36#include <gp_Circ2d.hxx>
37#include <gp_Torus.hxx>
38#include <gp_Cylinder.hxx>
39
40#include <Bnd_Box.hxx>
41
42#include <TColgp_HArray1OfPnt2d.hxx>
43#include <TColgp_SequenceOfPnt2d.hxx>
44#include <TColgp_Array1OfPnt.hxx>
45#include <TColgp_Array1OfPnt2d.hxx>
46
47#include <IntAna_QuadQuadGeo.hxx>
48
49#include <IntSurf_PntOn2S.hxx>
50#include <IntSurf_LineOn2S.hxx>
51#include <IntSurf_PntOn2S.hxx>
52#include <IntSurf_ListOfPntOn2S.hxx>
53#include <IntRes2d_Domain.hxx>
54#include <ProjLib_Plane.hxx>
55
56#include <IntPatch_GLine.hxx>
57#include <IntPatch_RLine.hxx>
58#include <IntPatch_WLine.hxx>
59#include <IntPatch_ALine.hxx>
60#include <IntPatch_ALineToWLine.hxx>
61
62#include <ElSLib.hxx>
63#include <ElCLib.hxx>
64
65#include <Extrema_ExtCC.hxx>
66#include <Extrema_POnCurv.hxx>
67#include <BndLib_AddSurface.hxx>
68
69#include <Adaptor3d_SurfacePtr.hxx>
70#include <Adaptor2d_HLine2d.hxx>
71
72#include <GeomAbs_SurfaceType.hxx>
73#include <GeomAbs_CurveType.hxx>
74
75#include <Geom_Surface.hxx>
76#include <Geom_Line.hxx>
77#include <Geom_Circle.hxx>
78#include <Geom_Ellipse.hxx>
79#include <Geom_Parabola.hxx>
80#include <Geom_Hyperbola.hxx>
81#include <Geom_TrimmedCurve.hxx>
82#include <Geom_BSplineCurve.hxx>
83#include <Geom_RectangularTrimmedSurface.hxx>
84#include <Geom_OffsetSurface.hxx>
85#include <Geom_Curve.hxx>
86#include <Geom_Conic.hxx>
87
88#include <Geom2d_TrimmedCurve.hxx>
89#include <Geom2d_BSplineCurve.hxx>
90#include <Geom2d_Line.hxx>
91#include <Geom2d_Curve.hxx>
92#include <Geom2d_Circle.hxx>
93
94#include <Geom2dAPI_InterCurveCurve.hxx>
95#include <Geom2dInt_GInter.hxx>
96#include <GeomAdaptor_Curve.hxx>
97#include <GeomAdaptor_HSurface.hxx>
98#include <GeomAdaptor_Surface.hxx>
99#include <GeomLib_CheckBSplineCurve.hxx>
100#include <GeomLib_Check2dBSplineCurve.hxx>
101
102#include <GeomInt_WLApprox.hxx>
103#include <GeomProjLib.hxx>
104#include <GeomAPI_ProjectPointOnSurf.hxx>
105#include <Geom2dAdaptor_Curve.hxx>
106#include <TopoDS.hxx>
107#include <TopoDS_Edge.hxx>
108#include <TopExp_Explorer.hxx>
109
110#include <BRep_Tool.hxx>
111#include <BRepTools.hxx>
112#include <BRepAdaptor_Surface.hxx>
113
114#include <BOPTColStd_Dump.hxx>
115
116#include <IntTools_Curve.hxx>
117#include <IntTools_Tools.hxx>
118#include <IntTools_Tools.hxx>
119#include <IntTools_TopolTool.hxx>
120#include <IntTools_PntOnFace.hxx>
121#include <IntTools_PntOn2Faces.hxx>
122#include <IntTools_Context.hxx>
123#include <IntSurf_ListIteratorOfListOfPntOn2S.hxx>
124
125static
126 void RefineVector(gp_Vec2d& aV2D);
127
128static
129 void DumpWLine(const Handle(IntPatch_WLine)& aWLine);
130//
131static
132 void TolR3d(const TopoDS_Face& ,
133 const TopoDS_Face& ,
134 Standard_Real& );
135static
136 Handle(Geom_Curve) MakeBSpline (const Handle(IntPatch_WLine)&,
137 const Standard_Integer,
138 const Standard_Integer);
139
140static
141 void Parameters(const Handle(GeomAdaptor_HSurface)&,
142 const Handle(GeomAdaptor_HSurface)&,
143 const gp_Pnt&,
144 Standard_Real&,
145 Standard_Real&,
146 Standard_Real&,
147 Standard_Real&);
148
149static
150 void BuildPCurves (Standard_Real f,Standard_Real l,Standard_Real& Tol,
151 const Handle (Geom_Surface)& S,
152 const Handle (Geom_Curve)& C,
153 Handle (Geom2d_Curve)& C2d);
154
155static
156 void CorrectSurfaceBoundaries(const TopoDS_Face& theFace,
157 const Standard_Real theTolerance,
158 Standard_Real& theumin,
159 Standard_Real& theumax,
160 Standard_Real& thevmin,
161 Standard_Real& thevmax);
162static
163 Standard_Boolean NotUseSurfacesForApprox
164 (const TopoDS_Face& aF1,
165 const TopoDS_Face& aF2,
166 const Handle(IntPatch_WLine)& WL,
167 const Standard_Integer ifprm,
168 const Standard_Integer ilprm);
169
170static
171 Handle(IntPatch_WLine) ComputePurgedWLine(const Handle(IntPatch_WLine)& theWLine);
172
173static
174 Standard_Real AdjustPeriodic(const Standard_Real theParameter,
175 const Standard_Real parmin,
176 const Standard_Real parmax,
177 const Standard_Real thePeriod,
178 Standard_Real& theOffset);
179
180static
181 Handle(Geom2d_BSplineCurve) MakeBSpline2d(const Handle(IntPatch_WLine)& theWLine,
182 const Standard_Integer ideb,
183 const Standard_Integer ifin,
184 const Standard_Boolean onFirst);
185
186static
187 Standard_Boolean DecompositionOfWLine(const Handle(IntPatch_WLine)& theWLine,
188 const Handle(GeomAdaptor_HSurface)& theSurface1,
189 const Handle(GeomAdaptor_HSurface)& theSurface2,
190 const TopoDS_Face& theFace1,
191 const TopoDS_Face& theFace2,
192 const IntTools_LineConstructor& theLConstructor,
193 const Standard_Boolean theAvoidLConstructor,
194 IntPatch_SequenceOfLine& theNewLines,
195 Standard_Real& theReachedTol3d,
196 const Handle(IntTools_Context)& );
197
198static
199 Standard_Boolean ParameterOutOfBoundary(const Standard_Real theParameter,
200 const Handle(Geom_Curve)& theCurve,
201 const TopoDS_Face& theFace1,
202 const TopoDS_Face& theFace2,
203 const Standard_Real theOtherParameter,
204 const Standard_Boolean bIncreasePar,
205 Standard_Real& theNewParameter,
206 const Handle(IntTools_Context)& );
207
208static
209 Standard_Boolean IsCurveValid(Handle(Geom2d_Curve)& thePCurve);
210
211static
212 Standard_Boolean IsPointOnBoundary(const Standard_Real theParameter,
213 const Standard_Real theFirstBoundary,
214 const Standard_Real theSecondBoundary,
215 const Standard_Real theResolution,
216 Standard_Boolean& IsOnFirstBoundary);
217static
218 Standard_Boolean FindPoint(const gp_Pnt2d& theFirstPoint,
219 const gp_Pnt2d& theLastPoint,
220 const Standard_Real theUmin,
221 const Standard_Real theUmax,
222 const Standard_Real theVmin,
223 const Standard_Real theVmax,
224 gp_Pnt2d& theNewPoint);
225
226
227static
228 Standard_Integer ComputeTangentZones( const Handle(GeomAdaptor_HSurface)& theSurface1,
229 const Handle(GeomAdaptor_HSurface)& theSurface2,
230 const TopoDS_Face& theFace1,
231 const TopoDS_Face& theFace2,
232 Handle(TColgp_HArray1OfPnt2d)& theResultOnS1,
233 Handle(TColgp_HArray1OfPnt2d)& theResultOnS2,
234 Handle(TColStd_HArray1OfReal)& theResultRadius,
235 const Handle(IntTools_Context)& );
236
237static
238 Standard_Boolean FindPoint(const gp_Pnt2d& theFirstPoint,
239 const gp_Pnt2d& theLastPoint,
240 const Standard_Real theUmin,
241 const Standard_Real theUmax,
242 const Standard_Real theVmin,
243 const Standard_Real theVmax,
244 const gp_Pnt2d& theTanZoneCenter,
245 const Standard_Real theZoneRadius,
246 Handle(GeomAdaptor_HSurface) theGASurface,
247 gp_Pnt2d& theNewPoint);
248
249static
250 Standard_Boolean IsInsideTanZone(const gp_Pnt2d& thePoint,
251 const gp_Pnt2d& theTanZoneCenter,
252 const Standard_Real theZoneRadius,
253 Handle(GeomAdaptor_HSurface) theGASurface);
254
255static
256 gp_Pnt2d AdjustByNeighbour(const gp_Pnt2d& theaNeighbourPoint,
257 const gp_Pnt2d& theOriginalPoint,
258 Handle(GeomAdaptor_HSurface) theGASurface);
259static
260 Standard_Boolean ApproxWithPCurves(const gp_Cylinder& theCyl,
261 const gp_Sphere& theSph);
262
263static void PerformPlanes(const Handle(GeomAdaptor_HSurface)& theS1,
264 const Handle(GeomAdaptor_HSurface)& theS2,
265 const Standard_Real TolAng,
266 const Standard_Real TolTang,
267 const Standard_Boolean theApprox1,
268 const Standard_Boolean theApprox2,
269 IntTools_SequenceOfCurves& theSeqOfCurve,
270 Standard_Boolean& theTangentFaces);
271
272static Standard_Boolean ClassifyLin2d(const Handle(GeomAdaptor_HSurface)& theS,
273 const gp_Lin2d& theLin2d,
274 const Standard_Real theTol,
275 Standard_Real& theP1,
276 Standard_Real& theP2);
277//
278static
279 void ApproxParameters(const Handle(GeomAdaptor_HSurface)& aHS1,
280 const Handle(GeomAdaptor_HSurface)& aHS2,
281 Standard_Integer& iDegMin,
282 Standard_Integer& iNbIter,
283 Standard_Integer& iDegMax);
284
285static
286 void Tolerances(const Handle(GeomAdaptor_HSurface)& aHS1,
287 const Handle(GeomAdaptor_HSurface)& aHS2,
288 Standard_Real& aTolArc,
289 Standard_Real& aTolTang,
290 Standard_Real& aUVMaxStep,
291 Standard_Real& aDeflection);
292
293static
294 Standard_Boolean SortTypes(const GeomAbs_SurfaceType aType1,
295 const GeomAbs_SurfaceType aType2);
296static
297 Standard_Integer IndexType(const GeomAbs_SurfaceType aType);
298
299//
300static
301 Standard_Real MaxSquareDistance (const Standard_Real aT,
302 const Handle(Geom_Curve)& aC3D,
303 const Handle(Geom2d_Curve)& aC2D1,
304 const Handle(Geom2d_Curve)& aC2D2,
305 const Handle(GeomAdaptor_HSurface) myHS1,
306 const Handle(GeomAdaptor_HSurface) myHS2,
307 const TopoDS_Face& aF1,
308 const TopoDS_Face& aF2,
309 const Handle(IntTools_Context)& aCtx);
310//
311static
312 Standard_Real FindMaxSquareDistance (const Standard_Real aA,
313 const Standard_Real aB,
314 const Standard_Real aEps,
315 const Handle(Geom_Curve)& aC3D,
316 const Handle(Geom2d_Curve)& aC2D1,
317 const Handle(Geom2d_Curve)& aC2D2,
318 const Handle(GeomAdaptor_HSurface)& myHS1,
319 const Handle(GeomAdaptor_HSurface)& myHS2,
320 const TopoDS_Face& aF1,
321 const TopoDS_Face& aF2,
322 const Handle(IntTools_Context)& aCtx);
323
324//=======================================================================
325//function :
326//purpose :
327//=======================================================================
328IntTools_FaceFace::IntTools_FaceFace()
329{
330 myTangentFaces=Standard_False;
331 //
332 myHS1 = new GeomAdaptor_HSurface ();
333 myHS2 = new GeomAdaptor_HSurface ();
334 myTolReached2d=0.;
335 myTolReached3d=0.;
336 SetParameters(Standard_True, Standard_True, Standard_True, 1.e-07);
337
338}
339//=======================================================================
340//function : SetContext
341//purpose :
342//=======================================================================
343void IntTools_FaceFace::SetContext(const Handle(IntTools_Context)& aContext)
344{
345 myContext=aContext;
346}
347//=======================================================================
348//function : Context
349//purpose :
350//=======================================================================
351const Handle(IntTools_Context)& IntTools_FaceFace::Context()const
352{
353 return myContext;
354}
355//=======================================================================
356//function : Face1
357//purpose :
358//=======================================================================
359const TopoDS_Face& IntTools_FaceFace::Face1() const
360{
361 return myFace1;
362}
363//=======================================================================
364//function : Face2
365//purpose :
366//=======================================================================
367const TopoDS_Face& IntTools_FaceFace::Face2() const
368{
369 return myFace2;
370}
371//=======================================================================
372//function : TangentFaces
373//purpose :
374//=======================================================================
375Standard_Boolean IntTools_FaceFace::TangentFaces() const
376{
377 return myTangentFaces;
378}
379//=======================================================================
380//function : Points
381//purpose :
382//=======================================================================
383const IntTools_SequenceOfPntOn2Faces& IntTools_FaceFace::Points() const
384{
385 return myPnts;
386}
387//=======================================================================
388//function : IsDone
389//purpose :
390//=======================================================================
391Standard_Boolean IntTools_FaceFace::IsDone() const
392{
393 return myIsDone;
394}
395//=======================================================================
396//function : TolReached3d
397//purpose :
398//=======================================================================
399Standard_Real IntTools_FaceFace::TolReached3d() const
400{
401 return myTolReached3d;
402}
403//=======================================================================
404//function : Lines
405//purpose : return lines of intersection
406//=======================================================================
407const IntTools_SequenceOfCurves& IntTools_FaceFace::Lines() const
408{
409 StdFail_NotDone_Raise_if
410 (!myIsDone,
411 "IntTools_FaceFace::Lines() => !myIntersector.IsDone()");
412 return mySeqOfCurve;
413}
414//=======================================================================
415//function : TolReached2d
416//purpose :
417//=======================================================================
418Standard_Real IntTools_FaceFace::TolReached2d() const
419{
420 return myTolReached2d;
421}
422// =======================================================================
423// function: SetParameters
424//
425// =======================================================================
426void IntTools_FaceFace::SetParameters(const Standard_Boolean ToApproxC3d,
427 const Standard_Boolean ToApproxC2dOnS1,
428 const Standard_Boolean ToApproxC2dOnS2,
429 const Standard_Real ApproximationTolerance)
430{
431 myApprox = ToApproxC3d;
432 myApprox1 = ToApproxC2dOnS1;
433 myApprox2 = ToApproxC2dOnS2;
434 myTolApprox = ApproximationTolerance;
435}
436//=======================================================================
437//function : SetList
438//purpose :
439//=======================================================================
440void IntTools_FaceFace::SetList(IntSurf_ListOfPntOn2S& aListOfPnts)
441{
442 myListOfPnts = aListOfPnts;
443}
444//=======================================================================
445//function : Perform
446//purpose : intersect surfaces of the faces
447//=======================================================================
448 void IntTools_FaceFace::Perform(const TopoDS_Face& aF1,
449 const TopoDS_Face& aF2)
450{
451 Standard_Boolean hasCone, RestrictLine, bTwoPlanes, bReverse;
452 Standard_Integer aNbLin, aNbPnts, i, NbLinPP;
453 Standard_Real TolArc, TolTang, Deflection, UVMaxStep;
454 Standard_Real umin, umax, vmin, vmax;
455 Standard_Real aTolF1, aTolF2;
456 GeomAbs_SurfaceType aType1, aType2;
457 Handle(Geom_Surface) S1, S2;
458 Handle(IntTools_TopolTool) dom1, dom2;
459 BRepAdaptor_Surface aBAS1, aBAS2;
460 //
461 if (myContext.IsNull()) {
462 myContext=new IntTools_Context;
463 }
464 //
465 mySeqOfCurve.Clear();
466 myTolReached2d=0.;
467 myTolReached3d=0.;
468 myIsDone = Standard_False;
469 myNbrestr=0;//?
470 hasCone = Standard_False;
471 bTwoPlanes = Standard_False;
472 //
473 myFace1=aF1;
474 myFace2=aF2;
475 //
476 aBAS1.Initialize(myFace1, Standard_False);
477 aBAS2.Initialize(myFace2, Standard_False);
478 aType1=aBAS1.GetType();
479 aType2=aBAS2.GetType();
480 //
481 bReverse=SortTypes(aType1, aType2);
482 if (bReverse) {
483 myFace1=aF2;
484 myFace2=aF1;
485 aType1=aBAS2.GetType();
486 aType2=aBAS1.GetType();
487 //
488 if (myListOfPnts.Extent()) {
489 Standard_Real aU1,aV1,aU2,aV2;
490 IntSurf_ListIteratorOfListOfPntOn2S aItP2S;
491 //
492 aItP2S.Initialize(myListOfPnts);
493 for (; aItP2S.More(); aItP2S.Next()){
494 IntSurf_PntOn2S& aP2S=aItP2S.Value();
495 aP2S.Parameters(aU1,aV1,aU2,aV2);
496 aP2S.SetValue(aU2,aV2,aU1,aV1);
497 }
498 }
499 }
500 //
501 S1=BRep_Tool::Surface(myFace1);
502 S2=BRep_Tool::Surface(myFace2);
503 //
504 aTolF1=BRep_Tool::Tolerance(myFace1);
505 aTolF2=BRep_Tool::Tolerance(myFace2);
506 //
507 TolArc= aTolF1 + aTolF2;
508 TolTang = TolArc;
509 //
510 NbLinPP = 0;
511 if(aType1==GeomAbs_Plane && aType2==GeomAbs_Plane){
512 bTwoPlanes = Standard_True;
513
514 BRepTools::UVBounds(myFace1, umin, umax, vmin, vmax);
515 myHS1->ChangeSurface().Load(S1, umin, umax, vmin, vmax);
516 //
517 BRepTools::UVBounds(myFace2, umin, umax, vmin, vmax);
518 myHS2->ChangeSurface().Load(S2, umin, umax, vmin, vmax);
519 Standard_Real TolAng = 1.e-8;
520 PerformPlanes(myHS1, myHS2, TolAng, TolTang, myApprox1, myApprox2,
521 mySeqOfCurve, myTangentFaces);
522
523 myIsDone = Standard_True;
524
525 if(!myTangentFaces) {
526 //
527 NbLinPP = mySeqOfCurve.Length();
528 if(NbLinPP) {
529 Standard_Real aTolFMax;
530 //
531 myTolReached3d = 1.e-7;
532 //
533 aTolFMax=Max(aTolF1, aTolF2);
534 //
535 if (aTolFMax>myTolReached3d) {
536 myTolReached3d=aTolFMax;
537 }
538 myTolReached2d = myTolReached3d;
539 //
540 if (bReverse) {
541 Handle(Geom2d_Curve) aC2D1, aC2D2;
542 //
543 aNbLin=mySeqOfCurve.Length();
544 for (i=1; i<=aNbLin; ++i) {
545 IntTools_Curve& aIC=mySeqOfCurve(i);
546 aC2D1=aIC.FirstCurve2d();
547 aC2D2=aIC.SecondCurve2d();
548 //
549 aIC.SetFirstCurve2d(aC2D2);
550 aIC.SetSecondCurve2d(aC2D1);
551 }
552 }
553 }
554 }
555 return;
556 }//if(aType1==GeomAbs_Plane && aType2==GeomAbs_Plane){
557 //
558 if (aType1==GeomAbs_Plane &&
559 (aType2==GeomAbs_Cylinder ||
560 aType2==GeomAbs_Cone ||
561 aType2==GeomAbs_Torus)) {
562 Standard_Real dU, dV;
563 // F1
564 BRepTools::UVBounds(myFace1, umin, umax, vmin, vmax);
565 dU=0.1*(umax-umin);
566 dV=0.1*(vmax-vmin);
567 umin=umin-dU;
568 umax=umax+dU;
569 vmin=vmin-dV;
570 vmax=vmax+dV;
571 myHS1->ChangeSurface().Load(S1, umin, umax, vmin, vmax);
572 // F2
573 BRepTools::UVBounds(myFace2, umin, umax, vmin, vmax);
574 CorrectSurfaceBoundaries(myFace2, (aTolF1 + aTolF2) * 2., umin, umax, vmin, vmax);
575 myHS2->ChangeSurface().Load(S2, umin, umax, vmin, vmax);
576 //
577 if( aType2==GeomAbs_Cone ) {
578 TolArc = 0.0001;
579 hasCone = Standard_True;
580 }
581 }
582 //
583 else if ((aType1==GeomAbs_Cylinder||
584 aType1==GeomAbs_Cone ||
585 aType1==GeomAbs_Torus) &&
586 aType2==GeomAbs_Plane) {
587 Standard_Real dU, dV;
588 //F1
589 BRepTools::UVBounds(myFace1, umin, umax, vmin, vmax);
590 CorrectSurfaceBoundaries(myFace1, (aTolF1 + aTolF2) * 2., umin, umax, vmin, vmax);
591 myHS1->ChangeSurface().Load(S1, umin, umax, vmin, vmax);
592 // F2
593 BRepTools::UVBounds(myFace2, umin, umax, vmin, vmax);
594 dU=0.1*(umax-umin);
595 dV=0.1*(vmax-vmin);
596 umin=umin-dU;
597 umax=umax+dU;
598 vmin=vmin-dV;
599 vmax=vmax+dV;
600 myHS2->ChangeSurface().Load(S2, umin, umax, vmin, vmax);
601 //
602 if( aType1==GeomAbs_Cone ) {
603 TolArc = 0.0001;
604 hasCone = Standard_True;
605 }
606 }
607
608 //
609 else {
610 BRepTools::UVBounds(myFace1, umin, umax, vmin, vmax);
611 //
612 CorrectSurfaceBoundaries(myFace1, (aTolF1 + aTolF2) * 2., umin, umax, vmin, vmax);
613 //
614 myHS1->ChangeSurface().Load(S1, umin, umax, vmin, vmax);
615 //
616 BRepTools::UVBounds(myFace2, umin, umax, vmin, vmax);
617 //
618 CorrectSurfaceBoundaries(myFace2, (aTolF1 + aTolF2) * 2., umin, umax, vmin, vmax);
619 //
620 myHS2->ChangeSurface().Load(S2, umin, umax, vmin, vmax);
621 }
622 //
623 dom1 = new IntTools_TopolTool(myHS1);
624 dom2 = new IntTools_TopolTool(myHS2);
625 //
626 myLConstruct.Load(dom1, dom2, myHS1, myHS2);
627 //
628 Deflection = (hasCone) ? 0.085 : 0.1;
629 UVMaxStep = 0.001;
630 //
631 Tolerances(myHS1, myHS2, TolArc, TolTang, UVMaxStep, Deflection);
632 //
633 myIntersector.SetTolerances(TolArc, TolTang, UVMaxStep, Deflection);
634 //
635 RestrictLine = Standard_False;
636 //
637 if((myHS1->IsUClosed() && !myHS1->IsUPeriodic()) ||
638 (myHS1->IsVClosed() && !myHS1->IsVPeriodic()) ||
639 (myHS2->IsUClosed() && !myHS2->IsUPeriodic()) ||
640 (myHS2->IsVClosed() && !myHS2->IsVPeriodic())) {
641 RestrictLine = Standard_True;
642 }
643 //
644 if(((aType1 != GeomAbs_BSplineSurface) &&
645 (aType1 != GeomAbs_BezierSurface) &&
646 (aType1 != GeomAbs_OtherSurface)) &&
647 ((aType2 != GeomAbs_BSplineSurface) &&
648 (aType2 != GeomAbs_BezierSurface) &&
649 (aType2 != GeomAbs_OtherSurface))) {
650 RestrictLine = Standard_True;
651 //
652 if ((aType1 == GeomAbs_Torus) ||
653 (aType2 == GeomAbs_Torus) ) {
654 myListOfPnts.Clear();
655 }
656 }
657 //
658 if(!RestrictLine) {
659 TopExp_Explorer aExp;
660 //
661 for(i = 0; (!RestrictLine) && (i < 2); i++) {
662 const TopoDS_Face& aF=(!i) ? myFace1 : myFace2;
663 aExp.Init(aF, TopAbs_EDGE);
664 for(; aExp.More(); aExp.Next()) {
665 const TopoDS_Edge& aE=TopoDS::Edge(aExp.Current());
666 //
667 if(BRep_Tool::Degenerated(aE)) {
668 RestrictLine = Standard_True;
669 break;
670 }
671 }
672 }
673 }
674 //
675 myIntersector.Perform(myHS1, dom1, myHS2, dom2,
676 TolArc, TolTang,
677 myListOfPnts, RestrictLine);
678 //
679 myIsDone = myIntersector.IsDone();
680 if (myIsDone) {
681 myTangentFaces=myIntersector.TangentFaces();
682 if (myTangentFaces) {
683 return;
684 }
685 //
686 if(RestrictLine) {
687 myListOfPnts.Clear(); // to use LineConstructor
688 }
689 //
690 aNbLin = myIntersector.NbLines();
691 for (i=1; i<=aNbLin; ++i) {
692 MakeCurve(i, dom1, dom2);
693 }
694 //
695 ComputeTolReached3d();
696 //
697 if (bReverse) {
698 Handle(Geom2d_Curve) aC2D1, aC2D2;
699 //
700 aNbLin=mySeqOfCurve.Length();
701 for (i=1; i<=aNbLin; ++i) {
702 IntTools_Curve& aIC=mySeqOfCurve(i);
703 aC2D1=aIC.FirstCurve2d();
704 aC2D2=aIC.SecondCurve2d();
705 //
706 aIC.SetFirstCurve2d(aC2D2);
707 aIC.SetSecondCurve2d(aC2D1);
708 }
709 }
710 //
711 // Points
712 Standard_Real U1,V1,U2,V2;
713 IntTools_PntOnFace aPntOnF1, aPntOnF2;
714 IntTools_PntOn2Faces aPntOn2Faces;
715 //
716 aNbPnts=myIntersector.NbPnts();
717 for (i=1; i<=aNbPnts; ++i) {
718 const IntSurf_PntOn2S& aISPnt=myIntersector.Point(i).PntOn2S();
719 const gp_Pnt& aPnt=aISPnt.Value();
720 aISPnt.Parameters(U1,V1,U2,V2);
721 aPntOnF1.Init(myFace1, aPnt, U1, V1);
722 aPntOnF2.Init(myFace2, aPnt, U2, V2);
723 //
724 if (!bReverse) {
725 aPntOn2Faces.SetP1(aPntOnF1);
726 aPntOn2Faces.SetP2(aPntOnF2);
727 }
728 else {
729 aPntOn2Faces.SetP2(aPntOnF1);
730 aPntOn2Faces.SetP1(aPntOnF2);
731 }
732 myPnts.Append(aPntOn2Faces);
733 }
734 //
735 }
736}
737//=======================================================================
738//function :ComputeTolReached3d
739//purpose :
740//=======================================================================
741 void IntTools_FaceFace::ComputeTolReached3d()
742{
743 Standard_Integer aNbLin;
744 GeomAbs_SurfaceType aType1, aType2;
745 //
746 aNbLin=myIntersector.NbLines();
747 if (!aNbLin) {
748 return;
749 }
750 //
751 aType1=myHS1->Surface().GetType();
752 aType2=myHS2->Surface().GetType();
753 //
754 if (aType1==GeomAbs_Cylinder && aType2==GeomAbs_Cylinder) {
755 if (aNbLin==2){
756 Handle(IntPatch_Line) aIL1, aIL2;
757 IntPatch_IType aTL1, aTL2;
758 //
759 aIL1=myIntersector.Line(1);
760 aIL2=myIntersector.Line(2);
761 aTL1=aIL1->ArcType();
762 aTL2=aIL2->ArcType();
763 if (aTL1==IntPatch_Lin && aTL2==IntPatch_Lin) {
764 Standard_Real aD, aDTresh, dTol;
765 gp_Lin aL1, aL2;
766 //
767 dTol=1.e-8;
768 aDTresh=1.5e-6;
769 //
770 aL1=Handle(IntPatch_GLine)::DownCast(aIL1)->Line();
771 aL2=Handle(IntPatch_GLine)::DownCast(aIL2)->Line();
772 aD=aL1.Distance(aL2);
773 aD=0.5*aD;
774 if (aD<aDTresh) {
775 myTolReached3d=aD+dTol;
776 }
777 return;
778 }
779 }
780 //ZZ
781 if (aNbLin) {// Check the distances
782 Standard_Integer i, aNbP, j ;
783 Standard_Real aT1, aT2, dT, aD2, aD2Max, aEps, aT11, aT12;
784 //
785 aD2Max=0.;
786 aNbP=10;
787 aNbLin=mySeqOfCurve.Length();
788 //
789 for (i=1; i<=aNbLin; ++i) {
790 const IntTools_Curve& aIC=mySeqOfCurve(i);
791 const Handle(Geom_Curve)& aC3D=aIC.Curve();
792 const Handle(Geom2d_Curve)& aC2D1=aIC.FirstCurve2d();
793 const Handle(Geom2d_Curve)& aC2D2=aIC.SecondCurve2d();
794 //
795 if (aC3D.IsNull()) {
796 continue;
797 }
798 const Handle(Geom_BSplineCurve)& aBC=
799 Handle(Geom_BSplineCurve)::DownCast(aC3D);
800 if (aBC.IsNull()) {
801 continue;
802 }
803 //
804 aT1=aBC->FirstParameter();
805 aT2=aBC->LastParameter();
806 //
807 aEps=0.01*(aT2-aT1);
808 dT=(aT2-aT1)/aNbP;
809 for (j=1; j<aNbP; ++j) {
810 aT11=aT1+j*dT;
811 aT12=aT11+dT;
812 aD2=FindMaxSquareDistance(aT11, aT12, aEps, aC3D, aC2D1, aC2D2,
813 myHS1, myHS2, myFace1, myFace2, myContext);
814 if (aD2>aD2Max) {
815 aD2Max=aD2;
816 }
817 }
818 }//for (i=1; i<=aNbLin; ++i) {
819 //
820 myTolReached3d=sqrt(aD2Max);
821 }// if (aNbLin)
822 }// if (aType1==GeomAbs_Cylinder && aType2==GeomAbs_Cylinder) {
823 //
824 //904/G3 f
825 else if (aType1==GeomAbs_Plane && aType2==GeomAbs_Plane) {
826 Standard_Real aTolF1, aTolF2, aTolFMax, aTolTresh;
827 //
828 aTolTresh=1.e-7;
829 //
830 aTolF1 = BRep_Tool::Tolerance(myFace1);
831 aTolF2 = BRep_Tool::Tolerance(myFace2);
832 aTolFMax=Max(aTolF1, aTolF2);
833 //
834 if (aTolFMax>aTolTresh) {
835 myTolReached3d=aTolFMax;
836 }
837 }//if (aType1==GeomAbs_Plane && aType2==GeomAbs_Plane) {
838 //t
839 //IFV Bug OCC20297
840 else if((aType1 == GeomAbs_Cylinder && aType2 == GeomAbs_Plane) ||
841 (aType2 == GeomAbs_Cylinder && aType1 == GeomAbs_Plane)) {
842 if(aNbLin == 1) {
843 const Handle(IntPatch_Line)& aIL1 = myIntersector.Line(1);
844 if(aIL1->ArcType() == IntPatch_Circle) {
845 gp_Circ aCir = Handle(IntPatch_GLine)::DownCast(aIL1)->Circle();
846 gp_XYZ aCirDir = aCir.Axis().Direction().XYZ();
847 gp_XYZ aPlDir;
848 gp_Pln aPln;
849 if(aType1 == GeomAbs_Plane) {
850 aPln = myHS1->Surface().Plane();
851 }
852 else {
853 aPln = myHS2->Surface().Plane();
854 }
855 aPlDir = aPln.Axis().Direction().XYZ();
856 Standard_Real cs = aCirDir*aPlDir;
857 if(cs < 0.) aPlDir.Reverse();
858 Standard_Real eps = 1.e-14;
859 if(!aPlDir.IsEqual(aCirDir, eps)) {
860 Standard_Integer aNbP = 11;
861 Standard_Real dt = 2.*M_PI / (aNbP - 1), t;
862 for(t = 0.; t < 2.*M_PI; t += dt) {
863 Standard_Real d = aPln.Distance(ElCLib::Value(t, aCir));
864 if(myTolReached3d < d) myTolReached3d = d;
865 }
866 myTolReached3d *= 1.1;
867 }
868 } //aIL1->ArcType() == IntPatch_Circle
869 } //aNbLin == 1
870 } // aType1 == GeomAbs_Cylinder && aType2 == GeomAbs_Plane)
871 //End IFV Bug OCC20297
872 //
873 else if ((aType1==GeomAbs_Plane && aType2==GeomAbs_Torus) ||
874 (aType2==GeomAbs_Plane && aType1==GeomAbs_Torus)) {
875 aNbLin=mySeqOfCurve.Length();
876 if (aNbLin!=1) {
877 return;
878 }
879 //
880 Standard_Integer i, aNbP;
881 Standard_Real aT, aT1, aT2, dT, aUT, aVT, aUP, aVP;
882 Standard_Real aDP, aDT, aDmax;
883 gp_Pln aPln;
884 gp_Torus aTorus;
885 gp_Pnt aP, aPP, aPT;
886 //
887 const IntTools_Curve& aIC=mySeqOfCurve(1);
888 const Handle(Geom_Curve)& aC3D=aIC.Curve();
889 const Handle(Geom_BSplineCurve)& aBS=
890 Handle(Geom_BSplineCurve)::DownCast(aC3D);
891 if (aBS.IsNull()) {
892 return;
893 }
894 //
895 aT1=aBS->FirstParameter();
896 aT2=aBS->LastParameter();
897 //
898 aPln =(aType1==GeomAbs_Plane) ? myHS1->Plane() : myHS2->Plane();
899 aTorus=(aType1==GeomAbs_Plane) ? myHS2->Torus() : myHS1->Torus();
900 //
901 aDmax=-1.;
902 aNbP=11;
903 dT=(aT2-aT1)/(aNbP-1);
904 for (i=0; i<aNbP; ++i) {
905 aT=aT1+i*dT;
906 if (i==aNbP-1) {
907 aT=aT2;
908 }
909 //
910 aC3D->D0(aT, aP);
911 //
912 ElSLib::Parameters(aPln, aP, aUP, aVP);
913 aPP=ElSLib::Value(aUP, aVP, aPln);
914 aDP=aP.SquareDistance(aPP);
915 if (aDP>aDmax) {
916 aDmax=aDP;
917 }
918 //
919 ElSLib::Parameters(aTorus, aP, aUT, aVT);
920 aPT=ElSLib::Value(aUT, aVT, aTorus);
921 aDT=aP.SquareDistance(aPT);
922 if (aDT>aDmax) {
923 aDmax=aDT;
924 }
925 }
926 //
927 if (aDmax > myTolReached3d*myTolReached3d) {
928 myTolReached3d=sqrt(aDmax);
929 myTolReached3d=1.1*myTolReached3d;
930 }
931 }// if ((aType1==GeomAbs_Plane && aType2==GeomAbs_Torus) ||
932 //
933 else if ((aType1==GeomAbs_SurfaceOfRevolution && aType2==GeomAbs_Cylinder) ||
934 (aType2==GeomAbs_SurfaceOfRevolution && aType1==GeomAbs_Cylinder)) {
935 Standard_Integer i, j, aNbP;
936 Standard_Real aT, aT1, aT2, dT, aD2max, aD2;
937 //
938 aNbLin=mySeqOfCurve.Length();
939 aD2max=0.;
940 aNbP=11;
941 //
942 for (i=1; i<=aNbLin; ++i) {
943 const IntTools_Curve& aIC=mySeqOfCurve(i);
944 const Handle(Geom_Curve)& aC3D=aIC.Curve();
945 const Handle(Geom2d_Curve)& aC2D1=aIC.FirstCurve2d();
946 const Handle(Geom2d_Curve)& aC2D2=aIC.SecondCurve2d();
947 //
948 if (aC3D.IsNull()) {
949 continue;
950 }
951 const Handle(Geom_BSplineCurve)& aBC=
952 Handle(Geom_BSplineCurve)::DownCast(aC3D);
953 if (aBC.IsNull()) {
954 return;
955 }
956 //
957 aT1=aBC->FirstParameter();
958 aT2=aBC->LastParameter();
959 //
960 dT=(aT2-aT1)/(aNbP-1);
961 for (j=0; j<aNbP; ++j) {
962 aT=aT1+j*dT;
963 if (j==aNbP-1) {
964 aT=aT2;
965 }
966 //
967 aD2=MaxSquareDistance(aT, aC3D, aC2D1, aC2D2,
968 myHS1, myHS2, myFace1, myFace2, myContext);
969 if (aD2>aD2max) {
970 aD2max=aD2;
971 }
972 }//for (j=0; j<aNbP; ++j) {
973
974 }//for (i=1; i<=aNbLin; ++i) {
975 //
976 aD2=myTolReached3d*myTolReached3d;
977 if (aD2max > aD2) {
978 myTolReached3d=sqrt(aD2max);
979 }
980 }//if((aType1==GeomAbs_SurfaceOfRevolution ...
981 //modified by NIZNHY-PKV Thu Aug 30 13:31:10 2012f
982 else if ((aType1==GeomAbs_Plane && aType2==GeomAbs_Sphere) ||
983 (aType2==GeomAbs_Plane && aType1==GeomAbs_Sphere)) {
984 Standard_Integer i, j, aNbP;
985 Standard_Real aT, aT1, aT2, dT, aD2max, aD2, aEps, aT11, aT12;
986 //
987 aNbLin=mySeqOfCurve.Length();
988 aD2max=0.;
989 aNbP=10;
990 //
991 for (i=1; i<=aNbLin; ++i) {
992 const IntTools_Curve& aIC=mySeqOfCurve(i);
993 const Handle(Geom_Curve)& aC3D=aIC.Curve();
994 const Handle(Geom2d_Curve)& aC2D1=aIC.FirstCurve2d();
995 const Handle(Geom2d_Curve)& aC2D2=aIC.SecondCurve2d();
996 //
997 const Handle(Geom2d_BSplineCurve)& aBC2D1=
998 Handle(Geom2d_BSplineCurve)::DownCast(aC2D1);
999 const Handle(Geom2d_BSplineCurve)& aBC2D2=
1000 Handle(Geom2d_BSplineCurve)::DownCast(aC2D2);
1001 //
1002 if (aBC2D1.IsNull() && aBC2D2.IsNull()) {
1003 return;
1004 }
1005 //
1006 if (!aBC2D1.IsNull()) {
1007 aT1=aBC2D1->FirstParameter();
1008 aT2=aBC2D1->LastParameter();
1009 }
1010 else {
1011 aT1=aBC2D2->FirstParameter();
1012 aT2=aBC2D2->LastParameter();
1013 }
1014 //
1015 aEps=0.01*(aT2-aT1);
1016 dT=(aT2-aT1)/(aNbP-1);
1017 for (j=0; j<aNbP; ++j) {
1018 aT=aT1+j*dT;
1019 aT11=aT1+j*dT;
1020 aT12=aT11+dT;
1021 if (j==aNbP-1) {
1022 aT12=aT2;
1023 }
1024 //
1025 aD2=FindMaxSquareDistance(aT11, aT12, aEps, aC3D, aC2D1, aC2D2,
1026 myHS1, myHS2, myFace1, myFace2, myContext);
1027 if (aD2>aD2max) {
1028 aD2max=aD2;
1029 }
1030 }//for (j=0; j<aNbP; ++j) {
1031
1032 }//for (i=1; i<=aNbLin; ++i) {
1033 //
1034 aD2=myTolReached3d*myTolReached3d;
1035 if (aD2max > aD2) {
1036 myTolReached3d=sqrt(aD2max);
1037 }
1038 }//else if ((aType1==GeomAbs_Plane && aType2==GeomAbs_Sphere) ...
1039 //modified by NIZNHY-PKV Thu Aug 30 13:31:12 2012t
1040}
1041//=======================================================================
1042//function : MakeCurve
1043//purpose :
1044//=======================================================================
1045 void IntTools_FaceFace::MakeCurve(const Standard_Integer Index,
1046 const Handle(Adaptor3d_TopolTool)& dom1,
1047 const Handle(Adaptor3d_TopolTool)& dom2)
1048{
1049 Standard_Boolean bDone, rejectSurface, reApprox, bAvoidLineConstructor;
1050 Standard_Boolean ok;
1051 Standard_Integer i, j, aNbParts;
1052 Standard_Real fprm, lprm;
1053 Standard_Real Tolpc;
1054 Handle(IntPatch_Line) L;
1055 IntPatch_IType typl;
1056 Handle(Geom_Curve) newc;
1057 //
1058 const Standard_Real TOLCHECK =0.0000001;
1059 const Standard_Real TOLANGCHECK=0.1;
1060 //
1061 rejectSurface = Standard_False;
1062 reApprox = Standard_False;
1063
1064 reapprox:;
1065
1066 Tolpc = myTolApprox;
1067 bAvoidLineConstructor = Standard_False;
1068 L = myIntersector.Line(Index);
1069 typl = L->ArcType();
1070 //
1071 if(typl==IntPatch_Walking) {
1072 Handle(IntPatch_Line) anewL;
1073 //
1074 const Handle(IntPatch_WLine)& aWLine=
1075 Handle(IntPatch_WLine)::DownCast(L);
1076 //DEBf
1077 //DumpWLine(aWLine);
1078 //DEBt
1079 anewL = ComputePurgedWLine(aWLine);
1080 if(anewL.IsNull()) {
1081 return;
1082 }
1083 L = anewL;
1084 //DEBf
1085 /*
1086 { const Handle(IntPatch_WLine)& aWLineX=
1087 Handle(IntPatch_WLine)::DownCast(L);
1088 DumpWLine(aWLineX);
1089 }
1090 */
1091 //DEBt
1092 //
1093 if(!myListOfPnts.IsEmpty()) {
1094 bAvoidLineConstructor = Standard_True;
1095 }
1096
1097 Standard_Integer nbp = aWLine->NbPnts();
1098 const IntSurf_PntOn2S& p1 = aWLine->Point(1);
1099 const IntSurf_PntOn2S& p2 = aWLine->Point(nbp);
1100
1101 const gp_Pnt& P1 = p1.Value();
1102 const gp_Pnt& P2 = p2.Value();
1103
1104 if(P1.SquareDistance(P2) < 1.e-14) {
1105 bAvoidLineConstructor = Standard_False;
1106 }
1107
1108 }
1109 //
1110 // Line Constructor
1111 if(!bAvoidLineConstructor) {
1112 myLConstruct.Perform(L);
1113 //
1114 bDone=myLConstruct.IsDone();
1115 aNbParts=myLConstruct.NbParts();
1116 if (!bDone|| !aNbParts) {
1117 return;
1118 }
1119 }
1120 // Do the Curve
1121
1122
1123 typl=L->ArcType();
1124 switch (typl) {
1125 //########################################
1126 // Line, Parabola, Hyperbola
1127 //########################################
1128 case IntPatch_Lin:
1129 case IntPatch_Parabola:
1130 case IntPatch_Hyperbola: {
1131 if (typl == IntPatch_Lin) {
1132 newc =
1133 new Geom_Line (Handle(IntPatch_GLine)::DownCast(L)->Line());
1134 }
1135
1136 else if (typl == IntPatch_Parabola) {
1137 newc =
1138 new Geom_Parabola(Handle(IntPatch_GLine)::DownCast(L)->Parabola());
1139 }
1140
1141 else if (typl == IntPatch_Hyperbola) {
1142 newc =
1143 new Geom_Hyperbola (Handle(IntPatch_GLine)::DownCast(L)->Hyperbola());
1144 }
1145 //
1146 // myTolReached3d
1147 if (typl == IntPatch_Lin) {
1148 TolR3d (myFace1, myFace2, myTolReached3d);
1149 }
1150 //
1151 aNbParts=myLConstruct.NbParts();
1152 for (i=1; i<=aNbParts; i++) {
1153 myLConstruct.Part(i, fprm, lprm);
1154
1155 if (!Precision::IsNegativeInfinite(fprm) &&
1156 !Precision::IsPositiveInfinite(lprm)) {
1157 //
1158 IntTools_Curve aCurve;
1159 //
1160 Handle(Geom_TrimmedCurve) aCT3D=new Geom_TrimmedCurve(newc, fprm, lprm);
1161 aCurve.SetCurve(aCT3D);
1162 if (typl == IntPatch_Parabola) {
1163 Standard_Real aTolF1, aTolF2, aTolBase;
1164
1165 aTolF1 = BRep_Tool::Tolerance(myFace1);
1166 aTolF2 = BRep_Tool::Tolerance(myFace2);
1167 aTolBase=aTolF1+aTolF2;
1168 myTolReached3d=IntTools_Tools::CurveTolerance(aCT3D, aTolBase);
1169 }
1170 //
1171 aCurve.SetCurve(new Geom_TrimmedCurve(newc, fprm, lprm));
1172 if(myApprox1) {
1173 Handle (Geom2d_Curve) C2d;
1174 BuildPCurves(fprm, lprm, Tolpc, myHS1->ChangeSurface().Surface(), newc, C2d);
1175 if(Tolpc>myTolReached2d || myTolReached2d==0.) {
1176 myTolReached2d=Tolpc;
1177 }
1178 //
1179 aCurve.SetFirstCurve2d(new Geom2d_TrimmedCurve(C2d,fprm,lprm));
1180 }
1181 else {
1182 Handle(Geom2d_BSplineCurve) H1;
1183 //
1184 aCurve.SetFirstCurve2d(H1);
1185 }
1186
1187 if(myApprox2) {
1188 Handle (Geom2d_Curve) C2d;
1189 BuildPCurves(fprm,lprm,Tolpc,myHS2->ChangeSurface().Surface(),newc,C2d);
1190 if(Tolpc>myTolReached2d || myTolReached2d==0.) {
1191 myTolReached2d=Tolpc;
1192 }
1193 //
1194 aCurve.SetSecondCurve2d(new Geom2d_TrimmedCurve(C2d,fprm,lprm));
1195 }
1196 else {
1197 Handle(Geom2d_BSplineCurve) H1;
1198 //
1199 aCurve.SetSecondCurve2d(H1);
1200 }
1201 mySeqOfCurve.Append(aCurve);
1202 } // end of if (!Precision::IsNegativeInfinite(fprm) && !Precision::IsPositiveInfinite(lprm))
1203
1204 else {
1205 // on regarde si on garde
1206 //
1207 Standard_Boolean bFNIt, bLPIt;
1208 Standard_Real aTestPrm, dT=100.;
1209
1210 bFNIt=Precision::IsNegativeInfinite(fprm);
1211 bLPIt=Precision::IsPositiveInfinite(lprm);
1212
1213 aTestPrm=0.;
1214
1215 if (bFNIt && !bLPIt) {
1216 aTestPrm=lprm-dT;
1217 }
1218 else if (!bFNIt && bLPIt) {
1219 aTestPrm=fprm+dT;
1220 }
1221
1222 gp_Pnt ptref(newc->Value(aTestPrm));
1223 //
1224
1225 Standard_Real u1, v1, u2, v2, Tol;
1226
1227 Tol = Precision::Confusion();
1228 Parameters(myHS1, myHS2, ptref, u1, v1, u2, v2);
1229 ok = (dom1->Classify(gp_Pnt2d(u1, v1), Tol) != TopAbs_OUT);
1230 if(ok) {
1231 ok = (dom2->Classify(gp_Pnt2d(u2,v2),Tol) != TopAbs_OUT);
1232 }
1233 if (ok) {
1234 Handle(Geom2d_BSplineCurve) H1;
1235 mySeqOfCurve.Append(IntTools_Curve(newc, H1, H1));
1236 }
1237 }
1238 }// end of for (i=1; i<=myLConstruct.NbParts(); i++)
1239 }// case IntPatch_Lin: case IntPatch_Parabola: case IntPatch_Hyperbola:
1240 break;
1241
1242 //########################################
1243 // Circle and Ellipse
1244 //########################################
1245 case IntPatch_Circle:
1246 case IntPatch_Ellipse: {
1247
1248 if (typl == IntPatch_Circle) {
1249 newc = new Geom_Circle
1250 (Handle(IntPatch_GLine)::DownCast(L)->Circle());
1251 }
1252 else { //IntPatch_Ellipse
1253 newc = new Geom_Ellipse
1254 (Handle(IntPatch_GLine)::DownCast(L)->Ellipse());
1255 }
1256 //
1257 // myTolReached3d
1258 TolR3d (myFace1, myFace2, myTolReached3d);
1259 //
1260 aNbParts=myLConstruct.NbParts();
1261 //
1262 Standard_Real aPeriod, aNul;
1263 TColStd_SequenceOfReal aSeqFprm, aSeqLprm;
1264
1265 aNul=0.;
1266 aPeriod=M_PI+M_PI;
1267
1268 for (i=1; i<=aNbParts; i++) {
1269 myLConstruct.Part(i, fprm, lprm);
1270
1271 if (fprm < aNul && lprm > aNul) {
1272 // interval that goes through 0. is divided on two intervals;
1273 while (fprm<aNul || fprm>aPeriod) fprm=fprm+aPeriod;
1274 while (lprm<aNul || lprm>aPeriod) lprm=lprm+aPeriod;
1275 //
1276 if((aPeriod - fprm) > Tolpc) {
1277 aSeqFprm.Append(fprm);
1278 aSeqLprm.Append(aPeriod);
1279 }
1280 else {
1281 gp_Pnt P1 = newc->Value(fprm);
1282 gp_Pnt P2 = newc->Value(aPeriod);
1283 Standard_Real aTolDist = BRep_Tool::Tolerance(myFace1) + BRep_Tool::Tolerance(myFace2);
1284 aTolDist = (myTolReached3d > aTolDist) ? myTolReached3d : aTolDist;
1285
1286 if(P1.Distance(P2) > aTolDist) {
1287 Standard_Real anewpar = fprm;
1288
1289 if(ParameterOutOfBoundary(fprm, newc, myFace1, myFace2, lprm, Standard_False, anewpar, myContext)) {
1290 fprm = anewpar;
1291 }
1292 aSeqFprm.Append(fprm);
1293 aSeqLprm.Append(aPeriod);
1294 }
1295 }
1296
1297 //
1298 if((lprm - aNul) > Tolpc) {
1299 aSeqFprm.Append(aNul);
1300 aSeqLprm.Append(lprm);
1301 }
1302 else {
1303 gp_Pnt P1 = newc->Value(aNul);
1304 gp_Pnt P2 = newc->Value(lprm);
1305 Standard_Real aTolDist = BRep_Tool::Tolerance(myFace1) + BRep_Tool::Tolerance(myFace2);
1306 aTolDist = (myTolReached3d > aTolDist) ? myTolReached3d : aTolDist;
1307
1308 if(P1.Distance(P2) > aTolDist) {
1309 Standard_Real anewpar = lprm;
1310
1311 if(ParameterOutOfBoundary(lprm, newc, myFace1, myFace2, fprm, Standard_True, anewpar, myContext)) {
1312 lprm = anewpar;
1313 }
1314 aSeqFprm.Append(aNul);
1315 aSeqLprm.Append(lprm);
1316 }
1317 }
1318 }
1319 else {
1320 // usual interval
1321 aSeqFprm.Append(fprm);
1322 aSeqLprm.Append(lprm);
1323 }
1324 }
1325
1326 //
1327 aNbParts=aSeqFprm.Length();
1328 for (i=1; i<=aNbParts; i++) {
1329 fprm=aSeqFprm(i);
1330 lprm=aSeqLprm(i);
1331 //
1332 Standard_Real aRealEpsilon=RealEpsilon();
1333 if (Abs(fprm) > aRealEpsilon || Abs(lprm-2.*M_PI) > aRealEpsilon) {
1334 //==============================================
1335 ////
1336 IntTools_Curve aCurve;
1337 Handle(Geom_TrimmedCurve) aTC3D=new Geom_TrimmedCurve(newc,fprm,lprm);
1338 aCurve.SetCurve(aTC3D);
1339 fprm=aTC3D->FirstParameter();
1340 lprm=aTC3D->LastParameter ();
1341 ////
1342 if (typl == IntPatch_Circle || typl == IntPatch_Ellipse) {////
1343 if(myApprox1) {
1344 Handle (Geom2d_Curve) C2d;
1345 BuildPCurves(fprm,lprm,Tolpc,myHS1->ChangeSurface().Surface(),newc,C2d);
1346 if(Tolpc>myTolReached2d || myTolReached2d==0) {
1347 myTolReached2d=Tolpc;
1348 }
1349 //
1350 aCurve.SetFirstCurve2d(C2d);
1351 }
1352 else { ////
1353 Handle(Geom2d_BSplineCurve) H1;
1354 aCurve.SetFirstCurve2d(H1);
1355 }
1356
1357
1358 if(myApprox2) {
1359 Handle (Geom2d_Curve) C2d;
1360 BuildPCurves(fprm,lprm,Tolpc,myHS2->ChangeSurface().Surface(),newc,C2d);
1361 if(Tolpc>myTolReached2d || myTolReached2d==0) {
1362 myTolReached2d=Tolpc;
1363 }
1364 //
1365 aCurve.SetSecondCurve2d(C2d);
1366 }
1367 else {
1368 Handle(Geom2d_BSplineCurve) H1;
1369 aCurve.SetSecondCurve2d(H1);
1370 }
1371 }
1372
1373 else {
1374 Handle(Geom2d_BSplineCurve) H1;
1375 aCurve.SetFirstCurve2d(H1);
1376 aCurve.SetSecondCurve2d(H1);
1377 }
1378 mySeqOfCurve.Append(aCurve);
1379 //==============================================
1380 } //if (Abs(fprm) > RealEpsilon() || Abs(lprm-2.*M_PI) > RealEpsilon())
1381
1382 else {
1383 // on regarde si on garde
1384 //
1385 if (aNbParts==1) {
1386// if (Abs(fprm) < RealEpsilon() && Abs(lprm-2.*M_PI) < RealEpsilon()) {
1387 if (Abs(fprm) <= aRealEpsilon && Abs(lprm-2.*M_PI) <= aRealEpsilon) {
1388 IntTools_Curve aCurve;
1389 Handle(Geom_TrimmedCurve) aTC3D=new Geom_TrimmedCurve(newc,fprm,lprm);
1390 aCurve.SetCurve(aTC3D);
1391 fprm=aTC3D->FirstParameter();
1392 lprm=aTC3D->LastParameter ();
1393
1394 if(myApprox1) {
1395 Handle (Geom2d_Curve) C2d;
1396 BuildPCurves(fprm,lprm,Tolpc,myHS1->ChangeSurface().Surface(),newc,C2d);
1397 if(Tolpc>myTolReached2d || myTolReached2d==0) {
1398 myTolReached2d=Tolpc;
1399 }
1400 //
1401 aCurve.SetFirstCurve2d(C2d);
1402 }
1403 else { ////
1404 Handle(Geom2d_BSplineCurve) H1;
1405 aCurve.SetFirstCurve2d(H1);
1406 }
1407
1408 if(myApprox2) {
1409 Handle (Geom2d_Curve) C2d;
1410 BuildPCurves(fprm,lprm,Tolpc,myHS2->ChangeSurface().Surface(),newc,C2d);
1411 if(Tolpc>myTolReached2d || myTolReached2d==0) {
1412 myTolReached2d=Tolpc;
1413 }
1414 //
1415 aCurve.SetSecondCurve2d(C2d);
1416 }
1417 else {
1418 Handle(Geom2d_BSplineCurve) H1;
1419 aCurve.SetSecondCurve2d(H1);
1420 }
1421 mySeqOfCurve.Append(aCurve);
1422 break;
1423 }
1424 }
1425 //
1426 Standard_Real aTwoPIdiv17, u1, v1, u2, v2, Tol;
1427
1428 aTwoPIdiv17=2.*M_PI/17.;
1429
1430 for (j=0; j<=17; j++) {
1431 gp_Pnt ptref (newc->Value (j*aTwoPIdiv17));
1432 Tol = Precision::Confusion();
1433
1434 Parameters(myHS1, myHS2, ptref, u1, v1, u2, v2);
1435 ok = (dom1->Classify(gp_Pnt2d(u1,v1),Tol) != TopAbs_OUT);
1436 if(ok) {
1437 ok = (dom2->Classify(gp_Pnt2d(u2,v2),Tol) != TopAbs_OUT);
1438 }
1439 if (ok) {
1440 IntTools_Curve aCurve;
1441 aCurve.SetCurve(newc);
1442 //==============================================
1443 if (typl == IntPatch_Circle || typl == IntPatch_Ellipse) {
1444
1445 if(myApprox1) {
1446 Handle (Geom2d_Curve) C2d;
1447 BuildPCurves(fprm, lprm, Tolpc, myHS1->ChangeSurface().Surface(), newc, C2d);
1448 if(Tolpc>myTolReached2d || myTolReached2d==0) {
1449 myTolReached2d=Tolpc;
1450 }
1451 //
1452 aCurve.SetFirstCurve2d(C2d);
1453 }
1454 else {
1455 Handle(Geom2d_BSplineCurve) H1;
1456 aCurve.SetFirstCurve2d(H1);
1457 }
1458
1459 if(myApprox2) {
1460 Handle (Geom2d_Curve) C2d;
1461 BuildPCurves(fprm, lprm, Tolpc,myHS2->ChangeSurface().Surface(), newc, C2d);
1462 if(Tolpc>myTolReached2d || myTolReached2d==0) {
1463 myTolReached2d=Tolpc;
1464 }
1465 //
1466 aCurve.SetSecondCurve2d(C2d);
1467 }
1468
1469 else {
1470 Handle(Geom2d_BSplineCurve) H1;
1471 aCurve.SetSecondCurve2d(H1);
1472 }
1473 }// end of if (typl == IntPatch_Circle || typl == IntPatch_Ellipse)
1474
1475 else {
1476 Handle(Geom2d_BSplineCurve) H1;
1477 //
1478 aCurve.SetFirstCurve2d(H1);
1479 aCurve.SetSecondCurve2d(H1);
1480 }
1481 //==============================================
1482 //
1483 mySeqOfCurve.Append(aCurve);
1484 break;
1485
1486 }// end of if (ok) {
1487 }// end of for (Standard_Integer j=0; j<=17; j++)
1488 }// end of else { on regarde si on garde
1489 }// for (i=1; i<=myLConstruct.NbParts(); i++)
1490 }// IntPatch_Circle: IntPatch_Ellipse:
1491 break;
1492
1493 case IntPatch_Analytic: {
1494 IntSurf_Quadric quad1,quad2;
1495 GeomAbs_SurfaceType typs = myHS1->Surface().GetType();
1496
1497 switch (typs) {
1498 case GeomAbs_Plane:
1499 quad1.SetValue(myHS1->Surface().Plane());
1500 break;
1501 case GeomAbs_Cylinder:
1502 quad1.SetValue(myHS1->Surface().Cylinder());
1503 break;
1504 case GeomAbs_Cone:
1505 quad1.SetValue(myHS1->Surface().Cone());
1506 break;
1507 case GeomAbs_Sphere:
1508 quad1.SetValue(myHS1->Surface().Sphere());
1509 break;
1510 default:
1511 Standard_ConstructionError::Raise("GeomInt_IntSS::MakeCurve 1");
1512 }
1513
1514 typs = myHS2->Surface().GetType();
1515
1516 switch (typs) {
1517 case GeomAbs_Plane:
1518 quad2.SetValue(myHS2->Surface().Plane());
1519 break;
1520 case GeomAbs_Cylinder:
1521 quad2.SetValue(myHS2->Surface().Cylinder());
1522 break;
1523 case GeomAbs_Cone:
1524 quad2.SetValue(myHS2->Surface().Cone());
1525 break;
1526 case GeomAbs_Sphere:
1527 quad2.SetValue(myHS2->Surface().Sphere());
1528 break;
1529 default:
1530 Standard_ConstructionError::Raise("GeomInt_IntSS::MakeCurve 2");
1531 }
1532 //
1533 //=========
1534 IntPatch_ALineToWLine convert (quad1, quad2);
1535
1536 if (!myApprox) {
1537 aNbParts=myLConstruct.NbParts();
1538 for (i=1; i<=aNbParts; i++) {
1539 myLConstruct.Part(i, fprm, lprm);
1540 Handle(IntPatch_WLine) WL =
1541 convert.MakeWLine(Handle(IntPatch_ALine)::DownCast(L), fprm, lprm);
1542 //
1543 Handle(Geom2d_BSplineCurve) H1;
1544 Handle(Geom2d_BSplineCurve) H2;
1545
1546 if(myApprox1) {
1547 H1 = MakeBSpline2d(WL, 1, WL->NbPnts(), Standard_True);
1548 }
1549
1550 if(myApprox2) {
1551 H2 = MakeBSpline2d(WL, 1, WL->NbPnts(), Standard_False);
1552 }
1553 //
1554 mySeqOfCurve.Append(IntTools_Curve(MakeBSpline(WL,1,WL->NbPnts()), H1, H2));
1555 }
1556 } // if (!myApprox)
1557
1558 else { // myApprox=TRUE
1559 GeomInt_WLApprox theapp3d;
1560 //
1561 Standard_Real tol2d = myTolApprox;
1562 //
1563 theapp3d.SetParameters(myTolApprox, tol2d, 4, 8, 0, Standard_True);
1564
1565 aNbParts=myLConstruct.NbParts();
1566 for (i=1; i<=aNbParts; i++) {
1567 myLConstruct.Part(i, fprm, lprm);
1568 Handle(IntPatch_WLine) WL =
1569 convert.MakeWLine(Handle(IntPatch_ALine):: DownCast(L),fprm,lprm);
1570
1571 theapp3d.Perform(myHS1,myHS2,WL,Standard_True,myApprox1,myApprox2, 1, WL->NbPnts());
1572
1573 if (!theapp3d.IsDone()) {
1574 //
1575 Handle(Geom2d_BSplineCurve) H1;
1576 Handle(Geom2d_BSplineCurve) H2;
1577
1578 if(myApprox1) {
1579 H1 = MakeBSpline2d(WL, 1, WL->NbPnts(), Standard_True);
1580 }
1581
1582 if(myApprox2) {
1583 H2 = MakeBSpline2d(WL, 1, WL->NbPnts(), Standard_False);
1584 }
1585 //
1586 mySeqOfCurve.Append(IntTools_Curve(MakeBSpline(WL,1,WL->NbPnts()), H1, H2));
1587 }
1588
1589 else {
1590 if(myApprox1 || myApprox2) {
1591 if( theapp3d.TolReached2d()>myTolReached2d || myTolReached2d==0) {
1592 myTolReached2d = theapp3d.TolReached2d();
1593 }
1594 }
1595
1596 if( theapp3d.TolReached3d()>myTolReached3d || myTolReached3d==0) {
1597 myTolReached3d = theapp3d.TolReached3d();
1598 }
1599
1600 Standard_Integer aNbMultiCurves, nbpoles;
1601 aNbMultiCurves=theapp3d.NbMultiCurves();
1602 for (j=1; j<=aNbMultiCurves; j++) {
1603 const AppParCurves_MultiBSpCurve& mbspc = theapp3d.Value(j);
1604 nbpoles = mbspc.NbPoles();
1605
1606 TColgp_Array1OfPnt tpoles(1, nbpoles);
1607 mbspc.Curve(1, tpoles);
1608 Handle(Geom_BSplineCurve) BS=new Geom_BSplineCurve(tpoles,
1609 mbspc.Knots(),
1610 mbspc.Multiplicities(),
1611 mbspc.Degree());
1612
1613 GeomLib_CheckBSplineCurve Check(BS,TOLCHECK,TOLANGCHECK);
1614 Check.FixTangent(Standard_True,Standard_True);
1615 //
1616 IntTools_Curve aCurve;
1617 aCurve.SetCurve(BS);
1618
1619 if(myApprox1) {
1620 TColgp_Array1OfPnt2d tpoles2d(1,nbpoles);
1621 mbspc.Curve(2,tpoles2d);
1622 Handle(Geom2d_BSplineCurve) BS2=new Geom2d_BSplineCurve(tpoles2d,
1623 mbspc.Knots(),
1624 mbspc.Multiplicities(),
1625 mbspc.Degree());
1626
1627 GeomLib_Check2dBSplineCurve newCheck(BS2,TOLCHECK,TOLANGCHECK);
1628 newCheck.FixTangent(Standard_True,Standard_True);
1629 //
1630 aCurve.SetFirstCurve2d(BS2);
1631 }
1632 else {
1633 Handle(Geom2d_BSplineCurve) H1;
1634 aCurve.SetFirstCurve2d(H1);
1635 }
1636
1637 if(myApprox2) {
1638 TColgp_Array1OfPnt2d tpoles2d(1, nbpoles);
1639 Standard_Integer TwoOrThree;
1640 TwoOrThree=myApprox1 ? 3 : 2;
1641 mbspc.Curve(TwoOrThree, tpoles2d);
1642 Handle(Geom2d_BSplineCurve) BS2 =new Geom2d_BSplineCurve(tpoles2d,
1643 mbspc.Knots(),
1644 mbspc.Multiplicities(),
1645 mbspc.Degree());
1646
1647 GeomLib_Check2dBSplineCurve newCheck(BS2,TOLCHECK,TOLANGCHECK);
1648 newCheck.FixTangent(Standard_True,Standard_True);
1649 //
1650 aCurve.SetSecondCurve2d(BS2);
1651 }
1652 else {
1653 Handle(Geom2d_BSplineCurve) H2;
1654 aCurve.SetSecondCurve2d(H2);
1655 }
1656 //
1657 mySeqOfCurve.Append(aCurve);
1658
1659 }// for (j=1; j<=aNbMultiCurves; j++) {
1660 }// else from if (!theapp3d.IsDone())
1661 }// for (i=1; i<=aNbParts; i++) {
1662 }// else { // myApprox=TRUE
1663 }// case IntPatch_Analytic:
1664 break;
1665
1666 case IntPatch_Walking:{
1667 Handle(IntPatch_WLine) WL =
1668 Handle(IntPatch_WLine)::DownCast(L);
1669 //
1670 Standard_Integer ifprm, ilprm;
1671 //
1672 if (!myApprox) {
1673 aNbParts = 1;
1674 if(!bAvoidLineConstructor){
1675 aNbParts=myLConstruct.NbParts();
1676 }
1677 for (i=1; i<=aNbParts; ++i) {
1678 Handle(Geom2d_BSplineCurve) H1, H2;
1679 Handle(Geom_Curve) aBSp;
1680 //
1681 if(bAvoidLineConstructor) {
1682 ifprm = 1;
1683 ilprm = WL->NbPnts();
1684 }
1685 else {
1686 myLConstruct.Part(i, fprm, lprm);
1687 ifprm=(Standard_Integer)fprm;
1688 ilprm=(Standard_Integer)lprm;
1689 }
1690 //
1691 if(myApprox1) {
1692 H1 = MakeBSpline2d(WL, ifprm, ilprm, Standard_True);
1693 }
1694 //
1695 if(myApprox2) {
1696 H2 = MakeBSpline2d(WL, ifprm, ilprm, Standard_False);
1697 }
1698 //
1699 aBSp=MakeBSpline(WL, ifprm, ilprm);
1700 IntTools_Curve aIC(aBSp, H1, H2);
1701 mySeqOfCurve.Append(aIC);
1702 }// for (i=1; i<=aNbParts; ++i) {
1703 }// if (!myApprox) {
1704 //
1705 else { // X
1706 Standard_Boolean bIsDecomposited;
1707 Standard_Integer nbiter, aNbSeqOfL;
1708 Standard_Real tol2d;
1709 IntPatch_SequenceOfLine aSeqOfL;
1710 GeomInt_WLApprox theapp3d;
1711 Approx_ParametrizationType aParType = Approx_ChordLength;
1712 //
1713 Standard_Boolean anApprox1 = myApprox1;
1714 Standard_Boolean anApprox2 = myApprox2;
1715
1716 tol2d = myTolApprox;
1717
1718 GeomAbs_SurfaceType typs1, typs2;
1719 typs1 = myHS1->Surface().GetType();
1720 typs2 = myHS2->Surface().GetType();
1721 Standard_Boolean anWithPC = Standard_True;
1722
1723 if(typs1 == GeomAbs_Cylinder && typs2 == GeomAbs_Sphere) {
1724 anWithPC =
1725 ApproxWithPCurves(myHS1->Surface().Cylinder(), myHS2->Surface().Sphere());
1726 }
1727 else if (typs1 == GeomAbs_Sphere && typs2 == GeomAbs_Cylinder) {
1728 anWithPC =
1729 ApproxWithPCurves(myHS2->Surface().Cylinder(), myHS1->Surface().Sphere());
1730 }
1731 if(!anWithPC) {
1732 //aParType = Approx_Centripetal;
1733 myTolApprox = 1.e-5;
1734 anApprox1 = Standard_False;
1735 anApprox2 = Standard_False;
1736 //
1737 tol2d = myTolApprox;
1738 }
1739
1740 if(myHS1 == myHS2) {
1741 //
1742 theapp3d.SetParameters(myTolApprox, tol2d, 4, 8, 0, Standard_False, aParType);
1743 rejectSurface = Standard_True;
1744 }
1745 else {
1746 if(reApprox && !rejectSurface)
1747 theapp3d.SetParameters(myTolApprox, tol2d, 4, 8, 0, Standard_False, aParType);
1748 else {
1749 Standard_Integer iDegMax, iDegMin, iNbIter;
1750 //
1751 ApproxParameters(myHS1, myHS2, iDegMin, iDegMax, iNbIter);
1752 theapp3d.SetParameters(myTolApprox, tol2d, iDegMin, iDegMax, iNbIter, Standard_True, aParType);
1753 //
1754 }
1755 }
1756 //
1757 Standard_Real aReachedTol = Precision::Confusion();
1758 bIsDecomposited=DecompositionOfWLine(WL,
1759 myHS1,
1760 myHS2,
1761 myFace1,
1762 myFace2,
1763 myLConstruct,
1764 bAvoidLineConstructor,
1765 aSeqOfL,
1766 aReachedTol,
1767 myContext);
1768 if ( bIsDecomposited && ( myTolReached3d < aReachedTol ) )
1769 myTolReached3d = aReachedTol;
1770
1771 //
1772 aNbSeqOfL=aSeqOfL.Length();
1773 //
1774 if (bIsDecomposited) {
1775 nbiter=aNbSeqOfL;
1776 }
1777 else {
1778 nbiter=1;
1779 aNbParts=1;
1780 if (!bAvoidLineConstructor) {
1781 aNbParts=myLConstruct.NbParts();
1782 nbiter=aNbParts;
1783 }
1784 }
1785 //
1786 // nbiter=(bIsDecomposited) ? aSeqOfL.Length() :
1787 // ((bAvoidLineConstructor) ? 1 :aNbParts);
1788 //
1789 for(i = 1; i <= nbiter; ++i) {
1790 if(bIsDecomposited) {
1791 WL = Handle(IntPatch_WLine)::DownCast(aSeqOfL.Value(i));
1792 ifprm = 1;
1793 ilprm = WL->NbPnts();
1794 }
1795 else {
1796 if(bAvoidLineConstructor) {
1797 ifprm = 1;
1798 ilprm = WL->NbPnts();
1799 }
1800 else {
1801 myLConstruct.Part(i, fprm, lprm);
1802 ifprm = (Standard_Integer)fprm;
1803 ilprm = (Standard_Integer)lprm;
1804 }
1805 }
1806 //-- lbr :
1807 //-- Si une des surfaces est un plan , on approxime en 2d
1808 //-- sur cette surface et on remonte les points 2d en 3d.
1809 if(typs1 == GeomAbs_Plane) {
1810 theapp3d.Perform(myHS1, myHS2, WL, Standard_False,Standard_True, myApprox2,ifprm,ilprm);
1811 }
1812 else if(typs2 == GeomAbs_Plane) {
1813 theapp3d.Perform(myHS1,myHS2,WL,Standard_False,myApprox1,Standard_True,ifprm,ilprm);
1814 }
1815 else {
1816 //
1817 if (myHS1 != myHS2){
1818 if ((typs1==GeomAbs_BezierSurface || typs1==GeomAbs_BSplineSurface) &&
1819 (typs2==GeomAbs_BezierSurface || typs2==GeomAbs_BSplineSurface)) {
1820
1821 theapp3d.SetParameters(myTolApprox, tol2d, 4, 8, 0, Standard_True, aParType);
1822
1823 Standard_Boolean bUseSurfaces;
1824 bUseSurfaces=NotUseSurfacesForApprox(myFace1, myFace2, WL, ifprm, ilprm);
1825 if (bUseSurfaces) {
1826 // ######
1827 rejectSurface = Standard_True;
1828 // ######
1829 theapp3d.SetParameters(myTolApprox, tol2d, 4, 8, 0, Standard_False, aParType);
1830 }
1831 }
1832 }
1833 //
1834 theapp3d.Perform(myHS1,myHS2,WL,Standard_True,anApprox1,anApprox2,ifprm,ilprm);
1835 }
1836
1837 if (!theapp3d.IsDone()) {
1838 //
1839 Handle(Geom2d_BSplineCurve) H1;
1840 //
1841 Handle(Geom_Curve) aBSp=MakeBSpline(WL,ifprm, ilprm);
1842 Handle(Geom2d_BSplineCurve) H2;
1843
1844 if(myApprox1) {
1845 H1 = MakeBSpline2d(WL, ifprm, ilprm, Standard_True);
1846 }
1847
1848 if(myApprox2) {
1849 H2 = MakeBSpline2d(WL, ifprm, ilprm, Standard_False);
1850 }
1851 //
1852 IntTools_Curve aIC(aBSp, H1, H2);
1853 mySeqOfCurve.Append(aIC);
1854 }
1855
1856 else {
1857 if(myApprox1 || myApprox2 || (typs1==GeomAbs_Plane || typs2==GeomAbs_Plane)) {
1858 if( theapp3d.TolReached2d()>myTolReached2d || myTolReached2d==0.) {
1859 myTolReached2d = theapp3d.TolReached2d();
1860 }
1861 }
1862 if(typs1==GeomAbs_Plane || typs2==GeomAbs_Plane) {
1863 myTolReached3d = myTolReached2d;
1864 //
1865 if (typs1==GeomAbs_Torus || typs2==GeomAbs_Torus) {
1866 if (myTolReached3d<1.e-6) {
1867 myTolReached3d = theapp3d.TolReached3d();
1868 myTolReached3d=1.e-6;
1869 }
1870 }
1871 //
1872 }
1873 else if( theapp3d.TolReached3d()>myTolReached3d || myTolReached3d==0.) {
1874 myTolReached3d = theapp3d.TolReached3d();
1875 }
1876
1877 Standard_Integer aNbMultiCurves, nbpoles;
1878 aNbMultiCurves=theapp3d.NbMultiCurves();
1879 for (j=1; j<=aNbMultiCurves; j++) {
1880 if(typs1 == GeomAbs_Plane) {
1881 const AppParCurves_MultiBSpCurve& mbspc = theapp3d.Value(j);
1882 nbpoles = mbspc.NbPoles();
1883
1884 TColgp_Array1OfPnt2d tpoles2d(1,nbpoles);
1885 TColgp_Array1OfPnt tpoles(1,nbpoles);
1886
1887 mbspc.Curve(1,tpoles2d);
1888 const gp_Pln& Pln = myHS1->Surface().Plane();
1889 //
1890 Standard_Integer ik;
1891 for(ik = 1; ik<= nbpoles; ik++) {
1892 tpoles.SetValue(ik,
1893 ElSLib::Value(tpoles2d.Value(ik).X(),
1894 tpoles2d.Value(ik).Y(),
1895 Pln));
1896 }
1897 //
1898 Handle(Geom_BSplineCurve) BS =
1899 new Geom_BSplineCurve(tpoles,
1900 mbspc.Knots(),
1901 mbspc.Multiplicities(),
1902 mbspc.Degree());
1903 GeomLib_CheckBSplineCurve Check(BS,TOLCHECK,TOLANGCHECK);
1904 Check.FixTangent(Standard_True, Standard_True);
1905 //
1906 IntTools_Curve aCurve;
1907 aCurve.SetCurve(BS);
1908
1909 if(myApprox1) {
1910 Handle(Geom2d_BSplineCurve) BS1 =
1911 new Geom2d_BSplineCurve(tpoles2d,
1912 mbspc.Knots(),
1913 mbspc.Multiplicities(),
1914 mbspc.Degree());
1915 GeomLib_Check2dBSplineCurve Check1(BS1,TOLCHECK,TOLANGCHECK);
1916 Check1.FixTangent(Standard_True,Standard_True);
1917 //
1918 // ############################################
1919 if(!rejectSurface && !reApprox) {
1920 Standard_Boolean isValid = IsCurveValid(BS1);
1921 if(!isValid) {
1922 reApprox = Standard_True;
1923 goto reapprox;
1924 }
1925 }
1926 // ############################################
1927 aCurve.SetFirstCurve2d(BS1);
1928 }
1929 else {
1930 Handle(Geom2d_BSplineCurve) H1;
1931 aCurve.SetFirstCurve2d(H1);
1932 }
1933
1934 if(myApprox2) {
1935 mbspc.Curve(2, tpoles2d);
1936
1937 Handle(Geom2d_BSplineCurve) BS2 = new Geom2d_BSplineCurve(tpoles2d,
1938 mbspc.Knots(),
1939 mbspc.Multiplicities(),
1940 mbspc.Degree());
1941 GeomLib_Check2dBSplineCurve newCheck(BS2,TOLCHECK,TOLANGCHECK);
1942 newCheck.FixTangent(Standard_True,Standard_True);
1943
1944 // ###########################################
1945 if(!rejectSurface && !reApprox) {
1946 Standard_Boolean isValid = IsCurveValid(BS2);
1947 if(!isValid) {
1948 reApprox = Standard_True;
1949 goto reapprox;
1950 }
1951 }
1952 // ###########################################
1953 //
1954 aCurve.SetSecondCurve2d(BS2);
1955 }
1956 else {
1957 Handle(Geom2d_BSplineCurve) H2;
1958 //
1959 aCurve.SetSecondCurve2d(H2);
1960 }
1961 //
1962 mySeqOfCurve.Append(aCurve);
1963 }
1964
1965 else if(typs2 == GeomAbs_Plane) {
1966 const AppParCurves_MultiBSpCurve& mbspc = theapp3d.Value(j);
1967 nbpoles = mbspc.NbPoles();
1968
1969 TColgp_Array1OfPnt2d tpoles2d(1,nbpoles);
1970 TColgp_Array1OfPnt tpoles(1,nbpoles);
1971 mbspc.Curve((myApprox1==Standard_True)? 2 : 1,tpoles2d);
1972 const gp_Pln& Pln = myHS2->Surface().Plane();
1973 //
1974 Standard_Integer ik;
1975 for(ik = 1; ik<= nbpoles; ik++) {
1976 tpoles.SetValue(ik,
1977 ElSLib::Value(tpoles2d.Value(ik).X(),
1978 tpoles2d.Value(ik).Y(),
1979 Pln));
1980
1981 }
1982 //
1983 Handle(Geom_BSplineCurve) BS=new Geom_BSplineCurve(tpoles,
1984 mbspc.Knots(),
1985 mbspc.Multiplicities(),
1986 mbspc.Degree());
1987 GeomLib_CheckBSplineCurve Check(BS,TOLCHECK,TOLANGCHECK);
1988 Check.FixTangent(Standard_True,Standard_True);
1989 //
1990 IntTools_Curve aCurve;
1991 aCurve.SetCurve(BS);
1992
1993 if(myApprox2) {
1994 Handle(Geom2d_BSplineCurve) BS1=new Geom2d_BSplineCurve(tpoles2d,
1995 mbspc.Knots(),
1996 mbspc.Multiplicities(),
1997 mbspc.Degree());
1998 GeomLib_Check2dBSplineCurve Check1(BS1,TOLCHECK,TOLANGCHECK);
1999 Check1.FixTangent(Standard_True,Standard_True);
2000 //
2001 // ###########################################
2002 if(!rejectSurface && !reApprox) {
2003 Standard_Boolean isValid = IsCurveValid(BS1);
2004 if(!isValid) {
2005 reApprox = Standard_True;
2006 goto reapprox;
2007 }
2008 }
2009 // ###########################################
2010 aCurve.SetSecondCurve2d(BS1);
2011 }
2012 else {
2013 Handle(Geom2d_BSplineCurve) H2;
2014 aCurve.SetSecondCurve2d(H2);
2015 }
2016
2017 if(myApprox1) {
2018 mbspc.Curve(1,tpoles2d);
2019 Handle(Geom2d_BSplineCurve) BS2=new Geom2d_BSplineCurve(tpoles2d,
2020 mbspc.Knots(),
2021 mbspc.Multiplicities(),
2022 mbspc.Degree());
2023 GeomLib_Check2dBSplineCurve Check2(BS2,TOLCHECK,TOLANGCHECK);
2024 Check2.FixTangent(Standard_True,Standard_True);
2025 //
2026 // ###########################################
2027 if(!rejectSurface && !reApprox) {
2028 Standard_Boolean isValid = IsCurveValid(BS2);
2029 if(!isValid) {
2030 reApprox = Standard_True;
2031 goto reapprox;
2032 }
2033 }
2034 // ###########################################
2035 aCurve.SetFirstCurve2d(BS2);
2036 }
2037 else {
2038 Handle(Geom2d_BSplineCurve) H1;
2039 //
2040 aCurve.SetFirstCurve2d(H1);
2041 }
2042 //
2043 mySeqOfCurve.Append(aCurve);
2044 }
2045 else {
2046 const AppParCurves_MultiBSpCurve& mbspc = theapp3d.Value(j);
2047 nbpoles = mbspc.NbPoles();
2048 TColgp_Array1OfPnt tpoles(1,nbpoles);
2049 mbspc.Curve(1,tpoles);
2050 Handle(Geom_BSplineCurve) BS=new Geom_BSplineCurve(tpoles,
2051 mbspc.Knots(),
2052 mbspc.Multiplicities(),
2053 mbspc.Degree());
2054 GeomLib_CheckBSplineCurve Check(BS,TOLCHECK,TOLANGCHECK);
2055 Check.FixTangent(Standard_True,Standard_True);
2056 //
2057 IntTools_Curve aCurve;
2058 aCurve.SetCurve(BS);
2059
2060 if(myApprox1) {
2061 if(anApprox1) {
2062 TColgp_Array1OfPnt2d tpoles2d(1,nbpoles);
2063 mbspc.Curve(2,tpoles2d);
2064 Handle(Geom2d_BSplineCurve) BS1=new Geom2d_BSplineCurve(tpoles2d,
2065 mbspc.Knots(),
2066 mbspc.Multiplicities(),
2067 mbspc.Degree());
2068 GeomLib_Check2dBSplineCurve newCheck(BS1,TOLCHECK,TOLANGCHECK);
2069 newCheck.FixTangent(Standard_True,Standard_True);
2070 //
2071 aCurve.SetFirstCurve2d(BS1);
2072 }
2073 else {
2074 Handle(Geom2d_BSplineCurve) BS1;
2075 fprm = BS->FirstParameter();
2076 lprm = BS->LastParameter();
2077
2078 Handle(Geom2d_Curve) C2d;
2079 Standard_Real aTol = myTolApprox;
2080 BuildPCurves(fprm, lprm, aTol, myHS1->ChangeSurface().Surface(), BS, C2d);
2081 BS1 = Handle(Geom2d_BSplineCurve)::DownCast(C2d);
2082 aCurve.SetFirstCurve2d(BS1);
2083 }
2084
2085 }
2086 else {
2087 Handle(Geom2d_BSplineCurve) H1;
2088 //
2089 aCurve.SetFirstCurve2d(H1);
2090 }
2091 if(myApprox2) {
2092 if(anApprox2) {
2093 TColgp_Array1OfPnt2d tpoles2d(1,nbpoles);
2094 mbspc.Curve((myApprox1==Standard_True)? 3 : 2,tpoles2d);
2095 Handle(Geom2d_BSplineCurve) BS2=new Geom2d_BSplineCurve(tpoles2d,
2096 mbspc.Knots(),
2097 mbspc.Multiplicities(),
2098 mbspc.Degree());
2099 GeomLib_Check2dBSplineCurve newCheck(BS2,TOLCHECK,TOLANGCHECK);
2100 newCheck.FixTangent(Standard_True,Standard_True);
2101 //
2102 aCurve.SetSecondCurve2d(BS2);
2103 }
2104 else {
2105 Handle(Geom2d_BSplineCurve) BS2;
2106 fprm = BS->FirstParameter();
2107 lprm = BS->LastParameter();
2108
2109 Handle(Geom2d_Curve) C2d;
2110 Standard_Real aTol = myTolApprox;
2111 BuildPCurves(fprm, lprm, aTol, myHS2->ChangeSurface().Surface(), BS, C2d);
2112 BS2 = Handle(Geom2d_BSplineCurve)::DownCast(C2d);
2113 aCurve.SetSecondCurve2d(BS2);
2114 }
2115
2116 }
2117 else {
2118 Handle(Geom2d_BSplineCurve) H2;
2119 //
2120 aCurve.SetSecondCurve2d(H2);
2121 }
2122 //
2123 mySeqOfCurve.Append(aCurve);
2124 }
2125 }
2126 }
2127 }
2128 }// else { // X
2129 }// case IntPatch_Walking:{
2130 break;
2131
2132 case IntPatch_Restriction:
2133 break;
2134
2135 }
2136}
2137
2138//=======================================================================
2139//function : BuildPCurves
2140//purpose :
2141//=======================================================================
2142 void BuildPCurves (Standard_Real f,
2143 Standard_Real l,
2144 Standard_Real& Tol,
2145 const Handle (Geom_Surface)& S,
2146 const Handle (Geom_Curve)& C,
2147 Handle (Geom2d_Curve)& C2d)
2148{
2149
2150 Standard_Real umin,umax,vmin,vmax;
2151 //
2152
2153 if (C2d.IsNull()) {
2154
2155 // in class ProjLib_Function the range of parameters is shrank by 1.e-09
2156 if((l - f) > 2.e-09) {
2157 C2d = GeomProjLib::Curve2d(C,f,l,S,Tol);
2158 //
2159 if (C2d.IsNull()) {
2160 // proj. a circle that goes through the pole on a sphere to the sphere
2161 Tol=Tol+1.e-7;
2162 C2d = GeomProjLib::Curve2d(C,f,l,S,Tol);
2163 }
2164 }
2165 else {
2166 if((l - f) > Epsilon(Abs(f))) {
2167 GeomAPI_ProjectPointOnSurf aProjector1, aProjector2;
2168 gp_Pnt P1 = C->Value(f);
2169 gp_Pnt P2 = C->Value(l);
2170 aProjector1.Init(P1, S);
2171 aProjector2.Init(P2, S);
2172
2173 if(aProjector1.IsDone() && aProjector2.IsDone()) {
2174 Standard_Real U=0., V=0.;
2175 aProjector1.LowerDistanceParameters(U, V);
2176 gp_Pnt2d p1(U, V);
2177
2178 aProjector2.LowerDistanceParameters(U, V);
2179 gp_Pnt2d p2(U, V);
2180
2181 if(p1.Distance(p2) > gp::Resolution()) {
2182 TColgp_Array1OfPnt2d poles(1,2);
2183 TColStd_Array1OfReal knots(1,2);
2184 TColStd_Array1OfInteger mults(1,2);
2185 poles(1) = p1;
2186 poles(2) = p2;
2187 knots(1) = f;
2188 knots(2) = l;
2189 mults(1) = mults(2) = 2;
2190
2191 C2d = new Geom2d_BSplineCurve(poles,knots,mults,1);
2192
2193 // compute reached tolerance.begin
2194 gp_Pnt PMid = C->Value((f + l) * 0.5);
2195 aProjector1.Perform(PMid);
2196
2197 if(aProjector1.IsDone()) {
2198 aProjector1.LowerDistanceParameters(U, V);
2199 gp_Pnt2d pmidproj(U, V);
2200 gp_Pnt2d pmidcurve2d = C2d->Value((f + l) * 0.5);
2201 Standard_Real adist = pmidcurve2d.Distance(pmidproj);
2202 Tol = (adist > Tol) ? adist : Tol;
2203 }
2204 // compute reached tolerance.end
2205 }
2206 }
2207 }
2208 }
2209 //
2210 S->Bounds(umin, umax, vmin, vmax);
2211
2212 if (S->IsUPeriodic() && !C2d.IsNull()) {
2213 // Recadre dans le domaine UV de la face
2214 Standard_Real period, U0, du, aEps;
2215
2216 du =0.0;
2217 aEps=Precision::PConfusion();
2218 period = S->UPeriod();
2219 gp_Pnt2d Pf = C2d->Value(f);
2220 U0=Pf.X();
2221 //
2222 gp_Pnt2d Pl = C2d->Value(l);
2223
2224 U0 = Min(Pl.X(), U0);
2225// while(U0-umin<aEps) {
2226 while(U0-umin<-aEps) {
2227 U0+=period;
2228 du+=period;
2229 }
2230 //
2231 while(U0-umax>aEps) {
2232 U0-=period;
2233 du-=period;
2234 }
2235 if (du != 0) {
2236 gp_Vec2d T1(du,0.);
2237 C2d->Translate(T1);
2238 }
2239 }
2240 }
2241 if (C2d.IsNull()) {
2242 BOPTColStd_Dump::PrintMessage("BuildPCurves()=> Echec ProjLib\n");
2243 }
2244}
2245
2246//=======================================================================
2247//function : Parameters
2248//purpose :
2249//=======================================================================
2250 void Parameters(const Handle(GeomAdaptor_HSurface)& HS1,
2251 const Handle(GeomAdaptor_HSurface)& HS2,
2252 const gp_Pnt& Ptref,
2253 Standard_Real& U1,
2254 Standard_Real& V1,
2255 Standard_Real& U2,
2256 Standard_Real& V2)
2257{
2258
2259 IntSurf_Quadric quad1,quad2;
2260 GeomAbs_SurfaceType typs = HS1->Surface().GetType();
2261
2262 switch (typs) {
2263 case GeomAbs_Plane:
2264 quad1.SetValue(HS1->Surface().Plane());
2265 break;
2266 case GeomAbs_Cylinder:
2267 quad1.SetValue(HS1->Surface().Cylinder());
2268 break;
2269 case GeomAbs_Cone:
2270 quad1.SetValue(HS1->Surface().Cone());
2271 break;
2272 case GeomAbs_Sphere:
2273 quad1.SetValue(HS1->Surface().Sphere());
2274 break;
2275 default:
2276 Standard_ConstructionError::Raise("GeomInt_IntSS::MakeCurve");
2277 }
2278
2279 typs = HS2->Surface().GetType();
2280 switch (typs) {
2281 case GeomAbs_Plane:
2282 quad2.SetValue(HS2->Surface().Plane());
2283 break;
2284 case GeomAbs_Cylinder:
2285 quad2.SetValue(HS2->Surface().Cylinder());
2286 break;
2287 case GeomAbs_Cone:
2288 quad2.SetValue(HS2->Surface().Cone());
2289 break;
2290 case GeomAbs_Sphere:
2291 quad2.SetValue(HS2->Surface().Sphere());
2292 break;
2293 default:
2294 Standard_ConstructionError::Raise("GeomInt_IntSS::MakeCurve");
2295 }
2296
2297 quad1.Parameters(Ptref,U1,V1);
2298 quad2.Parameters(Ptref,U2,V2);
2299}
2300
2301//=======================================================================
2302//function : MakeBSpline
2303//purpose :
2304//=======================================================================
2305Handle(Geom_Curve) MakeBSpline (const Handle(IntPatch_WLine)& WL,
2306 const Standard_Integer ideb,
2307 const Standard_Integer ifin)
2308{
2309 Standard_Integer i,nbpnt = ifin-ideb+1;
2310 TColgp_Array1OfPnt poles(1,nbpnt);
2311 TColStd_Array1OfReal knots(1,nbpnt);
2312 TColStd_Array1OfInteger mults(1,nbpnt);
2313 Standard_Integer ipidebm1;
2314 for(i=1,ipidebm1=i+ideb-1; i<=nbpnt;ipidebm1++, i++) {
2315 poles(i) = WL->Point(ipidebm1).Value();
2316 mults(i) = 1;
2317 knots(i) = i-1;
2318 }
2319 mults(1) = mults(nbpnt) = 2;
2320 return
2321 new Geom_BSplineCurve(poles,knots,mults,1);
2322}
2323//
2324
2325//=======================================================================
2326//function : MakeBSpline2d
2327//purpose :
2328//=======================================================================
2329Handle(Geom2d_BSplineCurve) MakeBSpline2d(const Handle(IntPatch_WLine)& theWLine,
2330 const Standard_Integer ideb,
2331 const Standard_Integer ifin,
2332 const Standard_Boolean onFirst)
2333{
2334 Standard_Integer i, nbpnt = ifin-ideb+1;
2335 TColgp_Array1OfPnt2d poles(1,nbpnt);
2336 TColStd_Array1OfReal knots(1,nbpnt);
2337 TColStd_Array1OfInteger mults(1,nbpnt);
2338 Standard_Integer ipidebm1;
2339
2340 for(i = 1, ipidebm1 = i+ideb-1; i <= nbpnt; ipidebm1++, i++) {
2341 Standard_Real U, V;
2342 if(onFirst)
2343 theWLine->Point(ipidebm1).ParametersOnS1(U, V);
2344 else
2345 theWLine->Point(ipidebm1).ParametersOnS2(U, V);
2346 poles(i).SetCoord(U, V);
2347 mults(i) = 1;
2348 knots(i) = i-1;
2349 }
2350 mults(1) = mults(nbpnt) = 2;
2351
2352 return new Geom2d_BSplineCurve(poles,knots,mults,1);
2353}
2354//=======================================================================
2355//function : PrepareLines3D
2356//purpose :
2357//=======================================================================
2358 void IntTools_FaceFace::PrepareLines3D(const Standard_Boolean bToSplit)
2359{
2360 Standard_Integer i, aNbCurves;
2361 GeomAbs_SurfaceType aType1, aType2;
2362 IntTools_SequenceOfCurves aNewCvs;
2363 //
2364 // 1. Treatment closed curves
2365 aNbCurves=mySeqOfCurve.Length();
2366 for (i=1; i<=aNbCurves; ++i) {
2367 const IntTools_Curve& aIC=mySeqOfCurve(i);
2368 //
2369 if (bToSplit) {
2370 Standard_Integer j, aNbC;
2371 IntTools_SequenceOfCurves aSeqCvs;
2372 //
2373 aNbC=IntTools_Tools::SplitCurve(aIC, aSeqCvs);
2374 if (aNbC) {
2375 for (j=1; j<=aNbC; ++j) {
2376 const IntTools_Curve& aICNew=aSeqCvs(j);
2377 aNewCvs.Append(aICNew);
2378 }
2379 }
2380 else {
2381 aNewCvs.Append(aIC);
2382 }
2383 }
2384 else {
2385 aNewCvs.Append(aIC);
2386 }
2387 }
2388 //
2389 // 2. Plane\Cone intersection when we had 4 curves
2390 aType1=myHS1->GetType();
2391 aType2=myHS2->GetType();
2392 aNbCurves=aNewCvs.Length();
2393 //
2394 if ((aType1==GeomAbs_Plane && aType2==GeomAbs_Cone) ||
2395 (aType2==GeomAbs_Plane && aType1==GeomAbs_Cone)) {
2396 if (aNbCurves==4) {
2397 GeomAbs_CurveType aCType1;
2398 //
2399 aCType1=aNewCvs(1).Type();
2400 if (aCType1==GeomAbs_Line) {
2401 IntTools_SequenceOfCurves aSeqIn, aSeqOut;
2402 //
2403 for (i=1; i<=aNbCurves; ++i) {
2404 const IntTools_Curve& aIC=aNewCvs(i);
2405 aSeqIn.Append(aIC);
2406 }
2407 //
2408 IntTools_Tools::RejectLines(aSeqIn, aSeqOut);
2409 //
2410 aNewCvs.Clear();
2411 aNbCurves=aSeqOut.Length();
2412 for (i=1; i<=aNbCurves; ++i) {
2413 const IntTools_Curve& aIC=aSeqOut(i);
2414 aNewCvs.Append(aIC);
2415 }
2416 }
2417 }
2418 }// if ((aType1==GeomAbs_Plane && aType2==GeomAbs_Cone)...
2419 //
2420 // 3. Fill mySeqOfCurve
2421 mySeqOfCurve.Clear();
2422 aNbCurves=aNewCvs.Length();
2423 for (i=1; i<=aNbCurves; ++i) {
2424 const IntTools_Curve& aIC=aNewCvs(i);
2425 mySeqOfCurve.Append(aIC);
2426 }
2427}
2428//=======================================================================
2429//function : CorrectSurfaceBoundaries
2430//purpose :
2431//=======================================================================
2432 void CorrectSurfaceBoundaries(const TopoDS_Face& theFace,
2433 const Standard_Real theTolerance,
2434 Standard_Real& theumin,
2435 Standard_Real& theumax,
2436 Standard_Real& thevmin,
2437 Standard_Real& thevmax)
2438{
2439 Standard_Boolean enlarge, isuperiodic, isvperiodic;
2440 Standard_Real uinf, usup, vinf, vsup, delta;
2441 GeomAbs_SurfaceType aType;
2442 Handle(Geom_Surface) aSurface;
2443 //
2444 aSurface = BRep_Tool::Surface(theFace);
2445 aSurface->Bounds(uinf, usup, vinf, vsup);
2446 delta = theTolerance;
2447 enlarge = Standard_False;
2448 //
2449 GeomAdaptor_Surface anAdaptorSurface(aSurface);
2450 //
2451 if(aSurface->IsKind(STANDARD_TYPE(Geom_RectangularTrimmedSurface))) {
2452 Handle(Geom_Surface) aBasisSurface =
2453 (Handle(Geom_RectangularTrimmedSurface)::DownCast(aSurface))->BasisSurface();
2454
2455 if(aBasisSurface->IsKind(STANDARD_TYPE(Geom_RectangularTrimmedSurface)) ||
2456 aBasisSurface->IsKind(STANDARD_TYPE(Geom_OffsetSurface))) {
2457 return;
2458 }
2459 }
2460 //
2461 if(aSurface->IsKind(STANDARD_TYPE(Geom_OffsetSurface))) {
2462 Handle(Geom_Surface) aBasisSurface =
2463 (Handle(Geom_OffsetSurface)::DownCast(aSurface))->BasisSurface();
2464
2465 if(aBasisSurface->IsKind(STANDARD_TYPE(Geom_RectangularTrimmedSurface)) ||
2466 aBasisSurface->IsKind(STANDARD_TYPE(Geom_OffsetSurface))) {
2467 return;
2468 }
2469 }
2470 //
2471 isuperiodic = anAdaptorSurface.IsUPeriodic();
2472 isvperiodic = anAdaptorSurface.IsVPeriodic();
2473 //
2474 aType=anAdaptorSurface.GetType();
2475 if((aType==GeomAbs_BezierSurface) ||
2476 (aType==GeomAbs_BSplineSurface) ||
2477 (aType==GeomAbs_SurfaceOfExtrusion) ||
2478 (aType==GeomAbs_SurfaceOfRevolution)) {
2479 enlarge=Standard_True;
2480 }
2481 //
2482 if(!isuperiodic && enlarge) {
2483
2484 if((theumin - uinf) > delta )
2485 theumin -= delta;
2486 else {
2487 theumin = uinf;
2488 }
2489
2490 if((usup - theumax) > delta )
2491 theumax += delta;
2492 else
2493 theumax = usup;
2494 }
2495 //
2496 if(!isvperiodic && enlarge) {
2497 if((thevmin - vinf) > delta ) {
2498 thevmin -= delta;
2499 }
2500 else {
2501 thevmin = vinf;
2502 }
2503 if((vsup - thevmax) > delta ) {
2504 thevmax += delta;
2505 }
2506 else {
2507 thevmax = vsup;
2508 }
2509 }
2510 //
2511 {
2512 Standard_Integer aNbP;
2513 Standard_Real aXP, dXfact, aXmid, aX1, aX2, aTolPA;
2514 //
2515 aTolPA=Precision::Angular();
2516 // U
2517 if (isuperiodic) {
2518 aXP=anAdaptorSurface.UPeriod();
2519 dXfact=theumax-theumin;
2520 if (dXfact-aTolPA>aXP) {
2521 aXmid=0.5*(theumax+theumin);
2522 aNbP=RealToInt(aXmid/aXP);
2523 if (aXmid<0.) {
2524 aNbP=aNbP-1;
2525 }
2526 aX1=aNbP*aXP;
2527 if (theumin>aTolPA) {
2528 aX1=theumin+aNbP*aXP;
2529 }
2530 aX2=aX1+aXP;
2531 if (theumin<aX1) {
2532 theumin=aX1;
2533 }
2534 if (theumax>aX2) {
2535 theumax=aX2;
2536 }
2537 }
2538 }
2539 // V
2540 if (isvperiodic) {
2541 aXP=anAdaptorSurface.VPeriod();
2542 dXfact=thevmax-thevmin;
2543 if (dXfact-aTolPA>aXP) {
2544 aXmid=0.5*(thevmax+thevmin);
2545 aNbP=RealToInt(aXmid/aXP);
2546 if (aXmid<0.) {
2547 aNbP=aNbP-1;
2548 }
2549 aX1=aNbP*aXP;
2550 if (thevmin>aTolPA) {
2551 aX1=thevmin+aNbP*aXP;
2552 }
2553 aX2=aX1+aXP;
2554 if (thevmin<aX1) {
2555 thevmin=aX1;
2556 }
2557 if (thevmax>aX2) {
2558 thevmax=aX2;
2559 }
2560 }
2561 }
2562 }
2563 //
2564 if(isuperiodic || isvperiodic) {
2565 Standard_Boolean correct = Standard_False;
2566 Standard_Boolean correctU = Standard_False;
2567 Standard_Boolean correctV = Standard_False;
2568 Bnd_Box2d aBox;
2569 TopExp_Explorer anExp;
2570
2571 for(anExp.Init(theFace, TopAbs_EDGE); anExp.More(); anExp.Next()) {
2572 if(BRep_Tool::IsClosed(TopoDS::Edge(anExp.Current()), theFace)) {
2573 correct = Standard_True;
2574 Standard_Real f, l;
2575 TopoDS_Edge anEdge = TopoDS::Edge(anExp.Current());
2576
2577 for(Standard_Integer i = 0; i < 2; i++) {
2578 if(i==0) {
2579 anEdge.Orientation(TopAbs_FORWARD);
2580 }
2581 else {
2582 anEdge.Orientation(TopAbs_REVERSED);
2583 }
2584 Handle(Geom2d_Curve) aCurve = BRep_Tool::CurveOnSurface(anEdge, theFace, f, l);
2585
2586 if(aCurve.IsNull()) {
2587 correct = Standard_False;
2588 break;
2589 }
2590 Handle(Geom2d_Line) aLine = Handle(Geom2d_Line)::DownCast(aCurve);
2591
2592 if(aLine.IsNull()) {
2593 correct = Standard_False;
2594 break;
2595 }
2596 gp_Dir2d anUDir(1., 0.);
2597 gp_Dir2d aVDir(0., 1.);
2598 Standard_Real anAngularTolerance = Precision::Angular();
2599
2600 correctU = correctU || aLine->Position().Direction().IsParallel(aVDir, anAngularTolerance);
2601 correctV = correctV || aLine->Position().Direction().IsParallel(anUDir, anAngularTolerance);
2602
2603 gp_Pnt2d pp1 = aCurve->Value(f);
2604 aBox.Add(pp1);
2605 gp_Pnt2d pp2 = aCurve->Value(l);
2606 aBox.Add(pp2);
2607 }
2608 if(!correct)
2609 break;
2610 }
2611 }
2612
2613 if(correct) {
2614 Standard_Real umin, vmin, umax, vmax;
2615 aBox.Get(umin, vmin, umax, vmax);
2616
2617 if(isuperiodic && correctU) {
2618
2619 if(theumin < umin)
2620 theumin = umin;
2621
2622 if(theumax > umax) {
2623 theumax = umax;
2624 }
2625 }
2626 if(isvperiodic && correctV) {
2627
2628 if(thevmin < vmin)
2629 thevmin = vmin;
2630 if(thevmax > vmax)
2631 thevmax = vmax;
2632 }
2633 }
2634 }
2635}
2636//
2637//
2638// The block is dedicated to determine whether WLine [ifprm, ilprm]
2639// crosses the degenerated zone on each given surface or not.
2640// If Yes -> We will not use info about surfaces during approximation
2641// because inside degenerated zone of the surface the approx. algo.
2642// uses wrong values of normal, etc., and resulting curve will have
2643// oscillations that we would not like to have.
2644// PKV Tue Feb 12 2002
2645
2646
2647static
2648 Standard_Boolean IsDegeneratedZone(const gp_Pnt2d& aP2d,
2649 const Handle(Geom_Surface)& aS,
2650 const Standard_Integer iDir);
2651static
2652 Standard_Boolean IsPointInDegeneratedZone(const IntSurf_PntOn2S& aP2S,
2653 const TopoDS_Face& aF1,
2654 const TopoDS_Face& aF2);
2655//=======================================================================
2656//function : NotUseSurfacesForApprox
2657//purpose :
2658//=======================================================================
2659Standard_Boolean NotUseSurfacesForApprox(const TopoDS_Face& aF1,
2660 const TopoDS_Face& aF2,
2661 const Handle(IntPatch_WLine)& WL,
2662 const Standard_Integer ifprm,
2663 const Standard_Integer ilprm)
2664{
2665 Standard_Boolean bPInDZ;
2666
2667 Handle(IntSurf_LineOn2S) aLineOn2S=WL->Curve();
2668
2669 const IntSurf_PntOn2S& aP2Sfprm=aLineOn2S->Value(ifprm);
2670 bPInDZ=IsPointInDegeneratedZone(aP2Sfprm, aF1, aF2);
2671 if (bPInDZ) {
2672 return bPInDZ;
2673 }
2674
2675 const IntSurf_PntOn2S& aP2Slprm=aLineOn2S->Value(ilprm);
2676 bPInDZ=IsPointInDegeneratedZone(aP2Slprm, aF1, aF2);
2677
2678 return bPInDZ;
2679}
2680//=======================================================================
2681//function : IsPointInDegeneratedZone
2682//purpose :
2683//=======================================================================
2684Standard_Boolean IsPointInDegeneratedZone(const IntSurf_PntOn2S& aP2S,
2685 const TopoDS_Face& aF1,
2686 const TopoDS_Face& aF2)
2687
2688{
2689 Standard_Boolean bFlag=Standard_True;
2690 Standard_Real US11, US12, VS11, VS12, US21, US22, VS21, VS22;
2691 Standard_Real U1, V1, U2, V2, aDelta, aD;
2692 gp_Pnt2d aP2d;
2693
2694 Handle(Geom_Surface)aS1 = BRep_Tool::Surface(aF1);
2695 aS1->Bounds(US11, US12, VS11, VS12);
2696 GeomAdaptor_Surface aGAS1(aS1);
2697
2698 Handle(Geom_Surface)aS2 = BRep_Tool::Surface(aF2);
2699 aS1->Bounds(US21, US22, VS21, VS22);
2700 GeomAdaptor_Surface aGAS2(aS2);
2701 //
2702 //const gp_Pnt& aP=aP2S.Value();
2703 aP2S.Parameters(U1, V1, U2, V2);
2704 //
2705 aDelta=1.e-7;
2706 // Check on Surf 1
2707 aD=aGAS1.UResolution(aDelta);
2708 aP2d.SetCoord(U1, V1);
2709 if (fabs(U1-US11) < aD) {
2710 bFlag=IsDegeneratedZone(aP2d, aS1, 1);
2711 if (bFlag) {
2712 return bFlag;
2713 }
2714 }
2715 if (fabs(U1-US12) < aD) {
2716 bFlag=IsDegeneratedZone(aP2d, aS1, 1);
2717 if (bFlag) {
2718 return bFlag;
2719 }
2720 }
2721 aD=aGAS1.VResolution(aDelta);
2722 if (fabs(V1-VS11) < aDelta) {
2723 bFlag=IsDegeneratedZone(aP2d, aS1, 2);
2724 if (bFlag) {
2725 return bFlag;
2726 }
2727 }
2728 if (fabs(V1-VS12) < aDelta) {
2729 bFlag=IsDegeneratedZone(aP2d, aS1, 2);
2730 if (bFlag) {
2731 return bFlag;
2732 }
2733 }
2734 // Check on Surf 2
2735 aD=aGAS2.UResolution(aDelta);
2736 aP2d.SetCoord(U2, V2);
2737 if (fabs(U2-US21) < aDelta) {
2738 bFlag=IsDegeneratedZone(aP2d, aS2, 1);
2739 if (bFlag) {
2740 return bFlag;
2741 }
2742 }
2743 if (fabs(U2-US22) < aDelta) {
2744 bFlag=IsDegeneratedZone(aP2d, aS2, 1);
2745 if (bFlag) {
2746 return bFlag;
2747 }
2748 }
2749 aD=aGAS2.VResolution(aDelta);
2750 if (fabs(V2-VS21) < aDelta) {
2751 bFlag=IsDegeneratedZone(aP2d, aS2, 2);
2752 if (bFlag) {
2753 return bFlag;
2754 }
2755 }
2756 if (fabs(V2-VS22) < aDelta) {
2757 bFlag=IsDegeneratedZone(aP2d, aS2, 2);
2758 if (bFlag) {
2759 return bFlag;
2760 }
2761 }
2762 return !bFlag;
2763}
2764
2765//=======================================================================
2766//function : IsDegeneratedZone
2767//purpose :
2768//=======================================================================
2769Standard_Boolean IsDegeneratedZone(const gp_Pnt2d& aP2d,
2770 const Handle(Geom_Surface)& aS,
2771 const Standard_Integer iDir)
2772{
2773 Standard_Boolean bFlag=Standard_True;
2774 Standard_Real US1, US2, VS1, VS2, dY, dX, d1, d2, dD;
2775 Standard_Real aXm, aYm, aXb, aYb, aXe, aYe;
2776 aS->Bounds(US1, US2, VS1, VS2);
2777
2778 gp_Pnt aPm, aPb, aPe;
2779
2780 aXm=aP2d.X();
2781 aYm=aP2d.Y();
2782
2783 aS->D0(aXm, aYm, aPm);
2784
2785 dX=1.e-5;
2786 dY=1.e-5;
2787 dD=1.e-12;
2788
2789 if (iDir==1) {
2790 aXb=aXm;
2791 aXe=aXm;
2792 aYb=aYm-dY;
2793 if (aYb < VS1) {
2794 aYb=VS1;
2795 }
2796 aYe=aYm+dY;
2797 if (aYe > VS2) {
2798 aYe=VS2;
2799 }
2800 aS->D0(aXb, aYb, aPb);
2801 aS->D0(aXe, aYe, aPe);
2802
2803 d1=aPm.Distance(aPb);
2804 d2=aPm.Distance(aPe);
2805 if (d1 < dD && d2 < dD) {
2806 return bFlag;
2807 }
2808 return !bFlag;
2809 }
2810 //
2811 else if (iDir==2) {
2812 aYb=aYm;
2813 aYe=aYm;
2814 aXb=aXm-dX;
2815 if (aXb < US1) {
2816 aXb=US1;
2817 }
2818 aXe=aXm+dX;
2819 if (aXe > US2) {
2820 aXe=US2;
2821 }
2822 aS->D0(aXb, aYb, aPb);
2823 aS->D0(aXe, aYe, aPe);
2824
2825 d1=aPm.Distance(aPb);
2826 d2=aPm.Distance(aPe);
2827 if (d1 < dD && d2 < dD) {
2828 return bFlag;
2829 }
2830 return !bFlag;
2831 }
2832 return !bFlag;
2833}
2834
2835//=========================================================================
2836// static function : ComputePurgedWLine
2837// purpose : Removes equal points (leave one of equal points) from theWLine
2838// and recompute vertex parameters.
2839// Returns new WLine or null WLine if the number
2840// of the points is less than 2.
2841//=========================================================================
2842Handle(IntPatch_WLine) ComputePurgedWLine(const Handle(IntPatch_WLine)& theWLine) {
2843
2844 Standard_Integer i, k, v, nb, nbvtx;
2845 Handle(IntPatch_WLine) aResult;
2846 nbvtx = theWLine->NbVertex();
2847 nb = theWLine->NbPnts();
2848 if (nb==2) {
2849 const IntSurf_PntOn2S& p1 = theWLine->Point(1);
2850 const IntSurf_PntOn2S& p2 = theWLine->Point(2);
2851 if(p1.Value().IsEqual(p2.Value(), gp::Resolution())) {
2852 return aResult;
2853 }
2854 }
2855 //
2856 Handle(IntPatch_WLine) aLocalWLine;
2857 Handle(IntPatch_WLine) aTmpWLine = theWLine;
2858 Handle(IntSurf_LineOn2S) aLineOn2S = new IntSurf_LineOn2S();
2859 aLocalWLine = new IntPatch_WLine(aLineOn2S, Standard_False);
2860 for(i = 1; i <= nb; i++) {
2861 aLineOn2S->Add(theWLine->Point(i));
2862 }
2863
2864 for(v = 1; v <= nbvtx; v++) {
2865 aLocalWLine->AddVertex(theWLine->Vertex(v));
2866 }
2867
2868 for(i = 1; i <= aLineOn2S->NbPoints(); i++) {
2869 Standard_Integer aStartIndex = i + 1;
2870 Standard_Integer anEndIndex = i + 5;
2871 nb = aLineOn2S->NbPoints();
2872 anEndIndex = (anEndIndex > nb) ? nb : anEndIndex;
2873
2874 if((aStartIndex > nb) || (anEndIndex <= 1)) {
2875 continue;
2876 }
2877 k = aStartIndex;
2878
2879 while(k <= anEndIndex) {
2880
2881 if(i != k) {
2882 IntSurf_PntOn2S p1 = aLineOn2S->Value(i);
2883 IntSurf_PntOn2S p2 = aLineOn2S->Value(k);
2884
2885 if(p1.Value().IsEqual(p2.Value(), gp::Resolution())) {
2886 aTmpWLine = aLocalWLine;
2887 aLocalWLine = new IntPatch_WLine(aLineOn2S, Standard_False);
2888
2889 for(v = 1; v <= aTmpWLine->NbVertex(); v++) {
2890 IntPatch_Point aVertex = aTmpWLine->Vertex(v);
2891 Standard_Integer avertexindex = (Standard_Integer)aVertex.ParameterOnLine();
2892
2893 if(avertexindex >= k) {
2894 aVertex.SetParameter(aVertex.ParameterOnLine() - 1.);
2895 }
2896 aLocalWLine->AddVertex(aVertex);
2897 }
2898 aLineOn2S->RemovePoint(k);
2899 anEndIndex--;
2900 continue;
2901 }
2902 }
2903 k++;
2904 }
2905 }
2906
2907 if(aLineOn2S->NbPoints() > 1) {
2908 aResult = aLocalWLine;
2909 }
2910 return aResult;
2911}
2912
2913//=======================================================================
2914//function : TolR3d
2915//purpose :
2916//=======================================================================
2917void TolR3d(const TopoDS_Face& aF1,
2918 const TopoDS_Face& aF2,
2919 Standard_Real& myTolReached3d)
2920{
2921 Standard_Real aTolF1, aTolF2, aTolFMax, aTolTresh;
2922
2923 aTolTresh=2.999999e-3;
2924 aTolF1 = BRep_Tool::Tolerance(aF1);
2925 aTolF2 = BRep_Tool::Tolerance(aF2);
2926 aTolFMax=Max(aTolF1, aTolF2);
2927
2928 if (aTolFMax>aTolTresh) {
2929 myTolReached3d=aTolFMax;
2930 }
2931}
2932//=======================================================================
2933//function : AdjustPeriodic
2934//purpose :
2935//=======================================================================
2936Standard_Real AdjustPeriodic(const Standard_Real theParameter,
2937 const Standard_Real parmin,
2938 const Standard_Real parmax,
2939 const Standard_Real thePeriod,
2940 Standard_Real& theOffset)
2941{
2942 Standard_Real aresult;
2943 //
2944 theOffset = 0.;
2945 aresult = theParameter;
2946 while(aresult < parmin) {
2947 aresult += thePeriod;
2948 theOffset += thePeriod;
2949 }
2950
2951 while(aresult > parmax) {
2952 aresult -= thePeriod;
2953 theOffset -= thePeriod;
2954 }
2955 return aresult;
2956}
2957//=======================================================================
2958//function : IsPointOnBoundary
2959//purpose :
2960//=======================================================================
2961Standard_Boolean IsPointOnBoundary(const Standard_Real theParameter,
2962 const Standard_Real theFirstBoundary,
2963 const Standard_Real theSecondBoundary,
2964 const Standard_Real theResolution,
2965 Standard_Boolean& IsOnFirstBoundary)
2966{
2967 Standard_Boolean bRet;
2968 Standard_Integer i;
2969 Standard_Real adist;
2970 //
2971 bRet=Standard_False;
2972 for(i = 0; i < 2; ++i) {
2973 IsOnFirstBoundary = (i == 0);
2974 if (IsOnFirstBoundary) {
2975 adist = fabs(theParameter - theFirstBoundary);
2976 }
2977 else {
2978 adist = fabs(theParameter - theSecondBoundary);
2979 }
2980 if(adist < theResolution) {
2981 return !bRet;
2982 }
2983 }
2984 return bRet;
2985}
2986// ------------------------------------------------------------------------------------------------
2987// static function: FindPoint
2988// purpose:
2989// ------------------------------------------------------------------------------------------------
2990Standard_Boolean FindPoint(const gp_Pnt2d& theFirstPoint,
2991 const gp_Pnt2d& theLastPoint,
2992 const Standard_Real theUmin,
2993 const Standard_Real theUmax,
2994 const Standard_Real theVmin,
2995 const Standard_Real theVmax,
2996 gp_Pnt2d& theNewPoint) {
2997
2998 gp_Vec2d aVec(theFirstPoint, theLastPoint);
2999 Standard_Integer i = 0, j = 0;
3000
3001 for(i = 0; i < 4; i++) {
3002 gp_Vec2d anOtherVec;
3003 gp_Vec2d anOtherVecNormal;
3004 gp_Pnt2d aprojpoint = theLastPoint;
3005
3006 if((i % 2) == 0) {
3007 anOtherVec.SetX(0.);
3008 anOtherVec.SetY(1.);
3009 anOtherVecNormal.SetX(1.);
3010 anOtherVecNormal.SetY(0.);
3011
3012 if(i < 2)
3013 aprojpoint.SetX(theUmin);
3014 else
3015 aprojpoint.SetX(theUmax);
3016 }
3017 else {
3018 anOtherVec.SetX(1.);
3019 anOtherVec.SetY(0.);
3020 anOtherVecNormal.SetX(0.);
3021 anOtherVecNormal.SetY(1.);
3022
3023 if(i < 2)
3024 aprojpoint.SetY(theVmin);
3025 else
3026 aprojpoint.SetY(theVmax);
3027 }
3028 gp_Vec2d anormvec = aVec;
3029 anormvec.Normalize();
3030 RefineVector(anormvec);
3031 Standard_Real adot1 = anormvec.Dot(anOtherVecNormal);
3032
3033 if(fabs(adot1) < Precision::Angular())
3034 continue;
3035 Standard_Real adist = 0.;
3036 Standard_Boolean bIsOut = Standard_False;
3037
3038 if((i % 2) == 0) {
3039 adist = (i < 2) ? fabs(theLastPoint.X() - theUmin) : fabs(theLastPoint.X() - theUmax);
3040 bIsOut = (i < 2) ? (theLastPoint.X() < theUmin) : (theLastPoint.X() > theUmax);
3041 }
3042 else {
3043 adist = (i < 2) ? fabs(theLastPoint.Y() - theVmin) : fabs(theLastPoint.Y() - theVmax);
3044 bIsOut = (i < 2) ? (theLastPoint.Y() < theVmin) : (theLastPoint.Y() > theVmax);
3045 }
3046 Standard_Real anoffset = adist * anOtherVec.Dot(anormvec) / adot1;
3047
3048 for(j = 0; j < 2; j++) {
3049 anoffset = (j == 0) ? anoffset : -anoffset;
3050 gp_Pnt2d acurpoint(aprojpoint.XY() + (anOtherVec.XY()*anoffset));
3051 gp_Vec2d acurvec(theLastPoint, acurpoint);
3052 if ( bIsOut )
3053 acurvec.Reverse();
3054
3055 Standard_Real aDotX, anAngleX;
3056 //
3057 aDotX = aVec.Dot(acurvec);
3058 anAngleX = aVec.Angle(acurvec);
3059 //
3060 if(aDotX > 0. && fabs(anAngleX) < Precision::PConfusion()) {
3061 if((i % 2) == 0) {
3062 if((acurpoint.Y() >= theVmin) &&
3063 (acurpoint.Y() <= theVmax)) {
3064 theNewPoint = acurpoint;
3065 return Standard_True;
3066 }
3067 }
3068 else {
3069 if((acurpoint.X() >= theUmin) &&
3070 (acurpoint.X() <= theUmax)) {
3071 theNewPoint = acurpoint;
3072 return Standard_True;
3073 }
3074 }
3075 }
3076 }
3077 }
3078 return Standard_False;
3079}
3080
3081
3082// ------------------------------------------------------------------------------------------------
3083// static function: FindPoint
3084// purpose: Find point on the boundary of radial tangent zone
3085// ------------------------------------------------------------------------------------------------
3086Standard_Boolean FindPoint(const gp_Pnt2d& theFirstPoint,
3087 const gp_Pnt2d& theLastPoint,
3088 const Standard_Real theUmin,
3089 const Standard_Real theUmax,
3090 const Standard_Real theVmin,
3091 const Standard_Real theVmax,
3092 const gp_Pnt2d& theTanZoneCenter,
3093 const Standard_Real theZoneRadius,
3094 Handle(GeomAdaptor_HSurface) theGASurface,
3095 gp_Pnt2d& theNewPoint) {
3096 theNewPoint = theLastPoint;
3097
3098 if ( !IsInsideTanZone( theLastPoint, theTanZoneCenter, theZoneRadius, theGASurface) )
3099 return Standard_False;
3100
3101 Standard_Real aUResolution = theGASurface->UResolution( theZoneRadius );
3102 Standard_Real aVResolution = theGASurface->VResolution( theZoneRadius );
3103
3104 Standard_Real aRadius = ( aUResolution < aVResolution ) ? aUResolution : aVResolution;
3105 gp_Ax22d anAxis( theTanZoneCenter, gp_Dir2d(1, 0), gp_Dir2d(0, 1) );
3106 gp_Circ2d aCircle( anAxis, aRadius );
3107
3108 //
3109 gp_Vec2d aDir( theLastPoint.XY() - theFirstPoint.XY() );
3110 Standard_Real aLength = aDir.Magnitude();
3111 if ( aLength <= gp::Resolution() )
3112 return Standard_False;
3113 gp_Lin2d aLine( theFirstPoint, aDir );
3114
3115 //
3116 Handle(Geom2d_Line) aCLine = new Geom2d_Line( aLine );
3117 Handle(Geom2d_TrimmedCurve) aC1 = new Geom2d_TrimmedCurve( aCLine, 0, aLength );
3118 Handle(Geom2d_Circle) aC2 = new Geom2d_Circle( aCircle );
3119
3120 Standard_Real aTol = aRadius * 0.001;
3121 aTol = ( aTol < Precision::PConfusion() ) ? Precision::PConfusion() : aTol;
3122
3123 Geom2dAPI_InterCurveCurve anIntersector;
3124 anIntersector.Init( aC1, aC2, aTol );
3125
3126 if ( anIntersector.NbPoints() == 0 )
3127 return Standard_False;
3128
3129 Standard_Boolean aFound = Standard_False;
3130 Standard_Real aMinDist = aLength * aLength;
3131 Standard_Integer i = 0;
3132 for ( i = 1; i <= anIntersector.NbPoints(); i++ ) {
3133 gp_Pnt2d aPInt = anIntersector.Point( i );
3134 if ( aPInt.SquareDistance( theFirstPoint ) < aMinDist ) {
3135 if ( ( aPInt.X() >= theUmin ) && ( aPInt.X() <= theUmax ) &&
3136 ( aPInt.Y() >= theVmin ) && ( aPInt.Y() <= theVmax ) ) {
3137 theNewPoint = aPInt;
3138 aFound = Standard_True;
3139 }
3140 }
3141 }
3142
3143 return aFound;
3144}
3145
3146// ------------------------------------------------------------------------------------------------
3147// static function: IsInsideTanZone
3148// purpose: Check if point is inside a radial tangent zone
3149// ------------------------------------------------------------------------------------------------
3150Standard_Boolean IsInsideTanZone(const gp_Pnt2d& thePoint,
3151 const gp_Pnt2d& theTanZoneCenter,
3152 const Standard_Real theZoneRadius,
3153 Handle(GeomAdaptor_HSurface) theGASurface) {
3154
3155 Standard_Real aUResolution = theGASurface->UResolution( theZoneRadius );
3156 Standard_Real aVResolution = theGASurface->VResolution( theZoneRadius );
3157 Standard_Real aRadiusSQR = ( aUResolution < aVResolution ) ? aUResolution : aVResolution;
3158 aRadiusSQR *= aRadiusSQR;
3159 if ( thePoint.SquareDistance( theTanZoneCenter ) <= aRadiusSQR )
3160 return Standard_True;
3161 return Standard_False;
3162}
3163
3164// ------------------------------------------------------------------------------------------------
3165// static function: CheckTangentZonesExist
3166// purpose: Check if tangent zone exists
3167// ------------------------------------------------------------------------------------------------
3168Standard_Boolean CheckTangentZonesExist( const Handle(GeomAdaptor_HSurface)& theSurface1,
3169 const Handle(GeomAdaptor_HSurface)& theSurface2 )
3170{
3171 if ( ( theSurface1->GetType() != GeomAbs_Torus ) ||
3172 ( theSurface2->GetType() != GeomAbs_Torus ) )
3173 return Standard_False;
3174
3175 gp_Torus aTor1 = theSurface1->Torus();
3176 gp_Torus aTor2 = theSurface2->Torus();
3177
3178 if ( aTor1.Location().Distance( aTor2.Location() ) > Precision::Confusion() )
3179 return Standard_False;
3180
3181 if ( ( fabs( aTor1.MajorRadius() - aTor2.MajorRadius() ) > Precision::Confusion() ) ||
3182 ( fabs( aTor1.MinorRadius() - aTor2.MinorRadius() ) > Precision::Confusion() ) )
3183 return Standard_False;
3184
3185 if ( ( aTor1.MajorRadius() < aTor1.MinorRadius() ) ||
3186 ( aTor2.MajorRadius() < aTor2.MinorRadius() ) )
3187 return Standard_False;
3188 return Standard_True;
3189}
3190
3191// ------------------------------------------------------------------------------------------------
3192// static function: ComputeTangentZones
3193// purpose:
3194// ------------------------------------------------------------------------------------------------
3195Standard_Integer ComputeTangentZones( const Handle(GeomAdaptor_HSurface)& theSurface1,
3196 const Handle(GeomAdaptor_HSurface)& theSurface2,
3197 const TopoDS_Face& theFace1,
3198 const TopoDS_Face& theFace2,
3199 Handle(TColgp_HArray1OfPnt2d)& theResultOnS1,
3200 Handle(TColgp_HArray1OfPnt2d)& theResultOnS2,
3201 Handle(TColStd_HArray1OfReal)& theResultRadius,
3202 const Handle(IntTools_Context)& aContext)
3203{
3204 Standard_Integer aResult = 0;
3205 if ( !CheckTangentZonesExist( theSurface1, theSurface2 ) )
3206 return aResult;
3207
3208
3209 TColgp_SequenceOfPnt2d aSeqResultS1, aSeqResultS2;
3210 TColStd_SequenceOfReal aSeqResultRad;
3211
3212 gp_Torus aTor1 = theSurface1->Torus();
3213 gp_Torus aTor2 = theSurface2->Torus();
3214
3215 gp_Ax2 anax1( aTor1.Location(), aTor1.Axis().Direction() );
3216 gp_Ax2 anax2( aTor2.Location(), aTor2.Axis().Direction() );
3217 Standard_Integer j = 0;
3218
3219 for ( j = 0; j < 2; j++ ) {
3220 Standard_Real aCoef = ( j == 0 ) ? -1 : 1;
3221 Standard_Real aRadius1 = fabs(aTor1.MajorRadius() + aCoef * aTor1.MinorRadius());
3222 Standard_Real aRadius2 = fabs(aTor2.MajorRadius() + aCoef * aTor2.MinorRadius());
3223
3224 gp_Circ aCircle1( anax1, aRadius1 );
3225 gp_Circ aCircle2( anax2, aRadius2 );
3226
3227 // roughly compute radius of tangent zone for perpendicular case
3228 Standard_Real aCriteria = Precision::Confusion() * 0.5;
3229
3230 Standard_Real aT1 = aCriteria;
3231 Standard_Real aT2 = aCriteria;
3232 if ( j == 0 ) {
3233 // internal tangency
3234 Standard_Real aR = ( aRadius1 > aTor2.MinorRadius() ) ? aRadius1 : aTor2.MinorRadius();
3235 //aT1 = aCriteria * aCriteria + aR * aR - ( aR - aCriteria ) * ( aR - aCriteria );
3236 aT1 = 2. * aR * aCriteria;
3237 aT2 = aT1;
3238 }
3239 else {
3240 // external tangency
3241 Standard_Real aRb = ( aRadius1 > aTor2.MinorRadius() ) ? aRadius1 : aTor2.MinorRadius();
3242 Standard_Real aRm = ( aRadius1 < aTor2.MinorRadius() ) ? aRadius1 : aTor2.MinorRadius();
3243 Standard_Real aDelta = aRb - aCriteria;
3244 aDelta *= aDelta;
3245 aDelta -= aRm * aRm;
3246 aDelta /= 2. * (aRb - aRm);
3247 aDelta -= 0.5 * (aRb - aRm);
3248
3249 aT1 = 2. * aRm * (aRm - aDelta);
3250 aT2 = aT1;
3251 }
3252 aCriteria = ( aT1 > aT2) ? aT1 : aT2;
3253 if ( aCriteria > 0 )
3254 aCriteria = sqrt( aCriteria );
3255
3256 if ( aCriteria > 0.5 * aTor1.MinorRadius() ) {
3257 // too big zone -> drop to minimum
3258 aCriteria = Precision::Confusion();
3259 }
3260
3261 GeomAdaptor_Curve aC1( new Geom_Circle(aCircle1) );
3262 GeomAdaptor_Curve aC2( new Geom_Circle(aCircle2) );
3263 Extrema_ExtCC anExtrema(aC1, aC2, 0, 2. * M_PI, 0, 2. * M_PI,
3264 Precision::PConfusion(), Precision::PConfusion());
3265
3266 if ( anExtrema.IsDone() ) {
3267
3268 Standard_Integer i = 0;
3269 for ( i = 1; i <= anExtrema.NbExt(); i++ ) {
3270 if ( anExtrema.SquareDistance(i) > aCriteria * aCriteria )
3271 continue;
3272
3273 Extrema_POnCurv P1, P2;
3274 anExtrema.Points( i, P1, P2 );
3275
3276 Standard_Boolean bFoundResult = Standard_True;
3277 gp_Pnt2d pr1, pr2;
3278
3279 Standard_Integer surfit = 0;
3280 for ( surfit = 0; surfit < 2; surfit++ ) {
3281 GeomAPI_ProjectPointOnSurf& aProjector =
3282 (surfit == 0) ? aContext->ProjPS(theFace1) : aContext->ProjPS(theFace2);
3283
3284 gp_Pnt aP3d = (surfit == 0) ? P1.Value() : P2.Value();
3285 aProjector.Perform(aP3d);
3286
3287 if(!aProjector.IsDone())
3288 bFoundResult = Standard_False;
3289 else {
3290 if(aProjector.LowerDistance() > aCriteria) {
3291 bFoundResult = Standard_False;
3292 }
3293 else {
3294 Standard_Real foundU = 0, foundV = 0;
3295 aProjector.LowerDistanceParameters(foundU, foundV);
3296 if ( surfit == 0 )
3297 pr1 = gp_Pnt2d( foundU, foundV );
3298 else
3299 pr2 = gp_Pnt2d( foundU, foundV );
3300 }
3301 }
3302 }
3303 if ( bFoundResult ) {
3304 aSeqResultS1.Append( pr1 );
3305 aSeqResultS2.Append( pr2 );
3306 aSeqResultRad.Append( aCriteria );
3307
3308 // torus is u and v periodic
3309 const Standard_Real twoPI = M_PI + M_PI;
3310 Standard_Real arr1tmp[2] = {pr1.X(), pr1.Y()};
3311 Standard_Real arr2tmp[2] = {pr2.X(), pr2.Y()};
3312
3313 // iteration on period bounds
3314 for ( Standard_Integer k1 = 0; k1 < 2; k1++ ) {
3315 Standard_Real aBound = ( k1 == 0 ) ? 0 : twoPI;
3316 Standard_Real aShift = ( k1 == 0 ) ? twoPI : -twoPI;
3317
3318 // iteration on surfaces
3319 for ( Standard_Integer k2 = 0; k2 < 2; k2++ ) {
3320 Standard_Real* arr1 = ( k2 == 0 ) ? arr1tmp : arr2tmp;
3321 Standard_Real* arr2 = ( k2 != 0 ) ? arr1tmp : arr2tmp;
3322 TColgp_SequenceOfPnt2d& aSeqS1 = ( k2 == 0 ) ? aSeqResultS1 : aSeqResultS2;
3323 TColgp_SequenceOfPnt2d& aSeqS2 = ( k2 != 0 ) ? aSeqResultS1 : aSeqResultS2;
3324
3325 if (fabs(arr1[0] - aBound) < Precision::PConfusion()) {
3326 aSeqS1.Append( gp_Pnt2d( arr1[0] + aShift, arr1[1] ) );
3327 aSeqS2.Append( gp_Pnt2d( arr2[0], arr2[1] ) );
3328 aSeqResultRad.Append( aCriteria );
3329 }
3330 if (fabs(arr1[1] - aBound) < Precision::PConfusion()) {
3331 aSeqS1.Append( gp_Pnt2d( arr1[0], arr1[1] + aShift) );
3332 aSeqS2.Append( gp_Pnt2d( arr2[0], arr2[1] ) );
3333 aSeqResultRad.Append( aCriteria );
3334 }
3335 }
3336 } //
3337 }
3338 }
3339 }
3340 }
3341 aResult = aSeqResultRad.Length();
3342
3343 if ( aResult > 0 ) {
3344 theResultOnS1 = new TColgp_HArray1OfPnt2d( 1, aResult );
3345 theResultOnS2 = new TColgp_HArray1OfPnt2d( 1, aResult );
3346 theResultRadius = new TColStd_HArray1OfReal( 1, aResult );
3347
3348 for ( Standard_Integer i = 1 ; i <= aResult; i++ ) {
3349 theResultOnS1->SetValue( i, aSeqResultS1.Value(i) );
3350 theResultOnS2->SetValue( i, aSeqResultS2.Value(i) );
3351 theResultRadius->SetValue( i, aSeqResultRad.Value(i) );
3352 }
3353 }
3354 return aResult;
3355}
3356
3357// ------------------------------------------------------------------------------------------------
3358// static function: AdjustByNeighbour
3359// purpose:
3360// ------------------------------------------------------------------------------------------------
3361gp_Pnt2d AdjustByNeighbour(const gp_Pnt2d& theaNeighbourPoint,
3362 const gp_Pnt2d& theOriginalPoint,
3363 Handle(GeomAdaptor_HSurface) theGASurface) {
3364
3365 gp_Pnt2d ap1 = theaNeighbourPoint;
3366 gp_Pnt2d ap2 = theOriginalPoint;
3367
3368 if ( theGASurface->IsUPeriodic() ) {
3369 Standard_Real aPeriod = theGASurface->UPeriod();
3370 gp_Pnt2d aPTest = ap2;
3371 Standard_Real aSqDistMin = 1.e+100;
3372
3373 for ( Standard_Integer pIt = -1; pIt <= 1; pIt++) {
3374 aPTest.SetX( theOriginalPoint.X() + aPeriod * pIt );
3375 Standard_Real dd = ap1.SquareDistance( aPTest );
3376
3377 if ( dd < aSqDistMin ) {
3378 ap2 = aPTest;
3379 aSqDistMin = dd;
3380 }
3381 }
3382 }
3383 if ( theGASurface->IsVPeriodic() ) {
3384 Standard_Real aPeriod = theGASurface->VPeriod();
3385 gp_Pnt2d aPTest = ap2;
3386 Standard_Real aSqDistMin = 1.e+100;
3387
3388 for ( Standard_Integer pIt = -1; pIt <= 1; pIt++) {
3389 aPTest.SetY( theOriginalPoint.Y() + aPeriod * pIt );
3390 Standard_Real dd = ap1.SquareDistance( aPTest );
3391
3392 if ( dd < aSqDistMin ) {
3393 ap2 = aPTest;
3394 aSqDistMin = dd;
3395 }
3396 }
3397 }
3398 return ap2;
3399}
3400
3401// ------------------------------------------------------------------------------------------------
3402//function: DecompositionOfWLine
3403// purpose:
3404// ------------------------------------------------------------------------------------------------
3405Standard_Boolean DecompositionOfWLine(const Handle(IntPatch_WLine)& theWLine,
3406 const Handle(GeomAdaptor_HSurface)& theSurface1,
3407 const Handle(GeomAdaptor_HSurface)& theSurface2,
3408 const TopoDS_Face& theFace1,
3409 const TopoDS_Face& theFace2,
3410 const IntTools_LineConstructor& theLConstructor,
3411 const Standard_Boolean theAvoidLConstructor,
3412 IntPatch_SequenceOfLine& theNewLines,
3413 Standard_Real& theReachedTol3d,
3414 const Handle(IntTools_Context)& aContext)
3415{
3416
3417 Standard_Boolean bRet, bAvoidLineConstructor;
3418 Standard_Integer aNbPnts, aNbParts;
3419 //
3420 bRet=Standard_False;
3421 aNbPnts=theWLine->NbPnts();
3422 bAvoidLineConstructor=theAvoidLConstructor;
3423 //
3424 if(!aNbPnts){
3425 return bRet;
3426 }
3427 if (!bAvoidLineConstructor) {
3428 aNbParts=theLConstructor.NbParts();
3429 if (!aNbParts) {
3430 return bRet;
3431 }
3432 }
3433 //
3434 Standard_Boolean bIsPrevPointOnBoundary, bIsPointOnBoundary, bIsCurrentPointOnBoundary;
3435 Standard_Integer nblines, pit, i, j;
3436 Standard_Real aTol;
3437 TColStd_Array1OfListOfInteger anArrayOfLines(1, aNbPnts);
3438 TColStd_Array1OfInteger anArrayOfLineType(1, aNbPnts);
3439 TColStd_ListOfInteger aListOfPointIndex;
3440
3441 Handle(TColgp_HArray1OfPnt2d) aTanZoneS1;
3442 Handle(TColgp_HArray1OfPnt2d) aTanZoneS2;
3443 Handle(TColStd_HArray1OfReal) aTanZoneRadius;
3444 Standard_Integer aNbZone = ComputeTangentZones( theSurface1, theSurface2, theFace1, theFace2,
3445 aTanZoneS1, aTanZoneS2, aTanZoneRadius, aContext);
3446
3447 //
3448 nblines=0;
3449 aTol=Precision::Confusion();
3450 aTol=0.5*aTol;
3451 bIsPrevPointOnBoundary=Standard_False;
3452 bIsPointOnBoundary=Standard_False;
3453 //
3454 // 1. ...
3455 //
3456 // Points
3457 for(pit = 1; pit <= aNbPnts; ++pit) {
3458 Standard_Boolean bIsOnFirstBoundary, isperiodic;
3459 Standard_Real aResolution, aPeriod, alowerboundary, aupperboundary, U, V;
3460 Standard_Real aParameter, anoffset, anAdjustPar;
3461 Standard_Real umin, umax, vmin, vmax;
3462 //
3463 bIsCurrentPointOnBoundary = Standard_False;
3464 const IntSurf_PntOn2S& aPoint = theWLine->Point(pit);
3465 //
3466 // Surface
3467 for(i = 0; i < 2; ++i) {
3468 Handle(GeomAdaptor_HSurface) aGASurface = (!i) ? theSurface1 : theSurface2;
3469 aGASurface->ChangeSurface().Surface()->Bounds(umin, umax, vmin, vmax);
3470 if(!i) {
3471 aPoint.ParametersOnS1(U, V);
3472 }
3473 else {
3474 aPoint.ParametersOnS2(U, V);
3475 }
3476 // U, V
3477 for(j = 0; j < 2; j++) {
3478 isperiodic = (!j) ? aGASurface->IsUPeriodic() : aGASurface->IsVPeriodic();
3479 if(!isperiodic){
3480 continue;
3481 }
3482 //
3483 if (!j) {
3484 aResolution=aGASurface->UResolution(aTol);
3485 aPeriod=aGASurface->UPeriod();
3486 alowerboundary=umin;
3487 aupperboundary=umax;
3488 aParameter=U;
3489 }
3490 else {
3491 aResolution=aGASurface->VResolution(aTol);
3492 aPeriod=aGASurface->VPeriod();
3493 alowerboundary=vmin;
3494 aupperboundary=vmax;
3495 aParameter=V;
3496 }
3497
3498 anoffset = 0.;
3499 anAdjustPar = AdjustPeriodic(aParameter,
3500 alowerboundary,
3501 aupperboundary,
3502 aPeriod,
3503 anoffset);
3504 //
3505 bIsOnFirstBoundary = Standard_True;// ?
3506 bIsPointOnBoundary=
3507 IsPointOnBoundary(anAdjustPar,
3508 alowerboundary,
3509 aupperboundary,
3510 aResolution,
3511 bIsOnFirstBoundary);
3512 //
3513 if(bIsPointOnBoundary) {
3514 bIsCurrentPointOnBoundary = Standard_True;
3515 break;
3516 }
3517 else {
3518 // check if a point belong to a tangent zone. Begin
3519 Standard_Integer zIt = 0;
3520 for ( zIt = 1; zIt <= aNbZone; zIt++ ) {
3521 gp_Pnt2d aPZone = (i == 0) ? aTanZoneS1->Value(zIt) : aTanZoneS2->Value(zIt);
3522 Standard_Real aZoneRadius = aTanZoneRadius->Value(zIt);
3523
3524 if ( IsInsideTanZone(gp_Pnt2d( U, V ), aPZone, aZoneRadius, aGASurface ) ) {
3525 // set boundary flag to split the curve by a tangent zone
3526 bIsPointOnBoundary = Standard_True;
3527 bIsCurrentPointOnBoundary = Standard_True;
3528 if ( theReachedTol3d < aZoneRadius ) {
3529 theReachedTol3d = aZoneRadius;
3530 }
3531 break;
3532 }
3533 }
3534 }
3535 }//for(j = 0; j < 2; j++) {
3536
3537 if(bIsCurrentPointOnBoundary){
3538 break;
3539 }
3540 }//for(i = 0; i < 2; ++i) {
3541 //
3542 if((bIsCurrentPointOnBoundary != bIsPrevPointOnBoundary)) {
3543 if(!aListOfPointIndex.IsEmpty()) {
3544 nblines++;
3545 anArrayOfLines.SetValue(nblines, aListOfPointIndex);
3546 anArrayOfLineType.SetValue(nblines, bIsPrevPointOnBoundary);
3547 aListOfPointIndex.Clear();
3548 }
3549 bIsPrevPointOnBoundary = bIsCurrentPointOnBoundary;
3550 }
3551 aListOfPointIndex.Append(pit);
3552 } //for(pit = 1; pit <= aNbPnts; ++pit) {
3553 //
3554 if(!aListOfPointIndex.IsEmpty()) {
3555 nblines++;
3556 anArrayOfLines.SetValue(nblines, aListOfPointIndex);
3557 anArrayOfLineType.SetValue(nblines, bIsPrevPointOnBoundary);
3558 aListOfPointIndex.Clear();
3559 }
3560 //
3561 if(nblines<=1) {
3562 return bRet; //Standard_False;
3563 }
3564 //
3565 //
3566 // 2. Correct wlines.begin
3567 TColStd_Array1OfListOfInteger anArrayOfLineEnds(1, nblines);
3568 Handle(IntSurf_LineOn2S) aSeqOfPntOn2S = new IntSurf_LineOn2S();
3569 //
3570 for(i = 1; i <= nblines; i++) {
3571 if(anArrayOfLineType.Value(i) != 0) {
3572 continue;
3573 }
3574 const TColStd_ListOfInteger& aListOfIndex = anArrayOfLines.Value(i);
3575 if(aListOfIndex.Extent() < 2) {
3576 continue;
3577 }
3578 TColStd_ListOfInteger aListOfFLIndex;
3579
3580 for(j = 0; j < 2; j++) {
3581 Standard_Integer aneighbourindex = (j == 0) ? (i - 1) : (i + 1);
3582
3583 if((aneighbourindex < 1) || (aneighbourindex > nblines))
3584 continue;
3585
3586 if(anArrayOfLineType.Value(aneighbourindex) == 0)
3587 continue;
3588 const TColStd_ListOfInteger& aNeighbour = anArrayOfLines.Value(aneighbourindex);
3589 Standard_Integer anIndex = (j == 0) ? aNeighbour.Last() : aNeighbour.First();
3590 const IntSurf_PntOn2S& aPoint = theWLine->Point(anIndex);
3591
3592 IntSurf_PntOn2S aNewP = aPoint;
3593
3594 for(Standard_Integer surfit = 0; surfit < 2; surfit++) {
3595
3596 Handle(GeomAdaptor_HSurface) aGASurface = (surfit == 0) ? theSurface1 : theSurface2;
3597 Standard_Real umin=0., umax=0., vmin=0., vmax=0.;
3598 aGASurface->ChangeSurface().Surface()->Bounds(umin, umax, vmin, vmax);
3599 Standard_Real U=0., V=0.;
3600
3601 if(surfit == 0)
3602 aNewP.ParametersOnS1(U, V);
3603 else
3604 aNewP.ParametersOnS2(U, V);
3605 Standard_Integer nbboundaries = 0;
3606
3607 Standard_Boolean bIsNearBoundary = Standard_False;
3608 Standard_Integer aZoneIndex = 0;
3609 Standard_Integer bIsUBoundary = Standard_False; // use if nbboundaries == 1
3610 Standard_Integer bIsFirstBoundary = Standard_False; // use if nbboundaries == 1
3611
3612
3613 for(Standard_Integer parit = 0; parit < 2; parit++) {
3614 Standard_Boolean isperiodic = (parit == 0) ? aGASurface->IsUPeriodic() : aGASurface->IsVPeriodic();
3615
3616 Standard_Real aResolution = (parit == 0) ? aGASurface->UResolution(aTol) : aGASurface->VResolution(aTol);
3617 Standard_Real alowerboundary = (parit == 0) ? umin : vmin;
3618 Standard_Real aupperboundary = (parit == 0) ? umax : vmax;
3619
3620 Standard_Real aParameter = (parit == 0) ? U : V;
3621 Standard_Boolean bIsOnFirstBoundary = Standard_True;
3622
3623 if(!isperiodic) {
3624 bIsPointOnBoundary=
3625 IsPointOnBoundary(aParameter, alowerboundary, aupperboundary, aResolution, bIsOnFirstBoundary);
3626 if(bIsPointOnBoundary) {
3627 bIsUBoundary = (parit == 0);
3628 bIsFirstBoundary = bIsOnFirstBoundary;
3629 nbboundaries++;
3630 }
3631 }
3632 else {
3633 Standard_Real aPeriod = (parit == 0) ? aGASurface->UPeriod() : aGASurface->VPeriod();
3634 Standard_Real anoffset = 0.;
3635 Standard_Real anAdjustPar = AdjustPeriodic(aParameter, alowerboundary, aupperboundary, aPeriod, anoffset);
3636
3637 bIsPointOnBoundary=
3638 IsPointOnBoundary(anAdjustPar, alowerboundary, aupperboundary, aResolution, bIsOnFirstBoundary);
3639 if(bIsPointOnBoundary) {
3640 bIsUBoundary = (parit == 0);
3641 bIsFirstBoundary = bIsOnFirstBoundary;
3642 nbboundaries++;
3643 }
3644 else {
3645 //check neighbourhood of boundary
3646 Standard_Real anEpsilon = aResolution * 100.;
3647 Standard_Real aPart = ( aupperboundary - alowerboundary ) * 0.1;
3648 anEpsilon = ( anEpsilon > aPart ) ? aPart : anEpsilon;
3649
3650 bIsNearBoundary = IsPointOnBoundary(anAdjustPar, alowerboundary, aupperboundary,
3651 anEpsilon, bIsOnFirstBoundary);
3652
3653 }
3654 }
3655 }
3656
3657 // check if a point belong to a tangent zone. Begin
3658 for ( Standard_Integer zIt = 1; zIt <= aNbZone; zIt++ ) {
3659 gp_Pnt2d aPZone = (surfit == 0) ? aTanZoneS1->Value(zIt) : aTanZoneS2->Value(zIt);
3660 Standard_Real aZoneRadius = aTanZoneRadius->Value(zIt);
3661
3662 Standard_Integer aneighbourpointindex1 = (j == 0) ? aListOfIndex.First() : aListOfIndex.Last();
3663 const IntSurf_PntOn2S& aNeighbourPoint = theWLine->Point(aneighbourpointindex1);
3664 Standard_Real nU1, nV1;
3665
3666 if(surfit == 0)
3667 aNeighbourPoint.ParametersOnS1(nU1, nV1);
3668 else
3669 aNeighbourPoint.ParametersOnS2(nU1, nV1);
3670 gp_Pnt2d ap1(nU1, nV1);
3671 gp_Pnt2d ap2 = AdjustByNeighbour( ap1, gp_Pnt2d( U, V ), aGASurface );
3672
3673
3674 if ( IsInsideTanZone( ap2, aPZone, aZoneRadius, aGASurface ) ) {
3675 aZoneIndex = zIt;
3676 bIsNearBoundary = Standard_True;
3677 if ( theReachedTol3d < aZoneRadius ) {
3678 theReachedTol3d = aZoneRadius;
3679 }
3680 }
3681 }
3682 // check if a point belong to a tangent zone. End
3683 Standard_Boolean bComputeLineEnd = Standard_False;
3684
3685 if(nbboundaries == 2) {
3686 //xf
3687 bComputeLineEnd = Standard_True;
3688 //xt
3689 }
3690 else if(nbboundaries == 1) {
3691 Standard_Boolean isperiodic = (bIsUBoundary) ? aGASurface->IsUPeriodic() : aGASurface->IsVPeriodic();
3692
3693 if(isperiodic) {
3694 Standard_Real alowerboundary = (bIsUBoundary) ? umin : vmin;
3695 Standard_Real aupperboundary = (bIsUBoundary) ? umax : vmax;
3696 Standard_Real aPeriod = (bIsUBoundary) ? aGASurface->UPeriod() : aGASurface->VPeriod();
3697 Standard_Real aParameter = (bIsUBoundary) ? U : V;
3698 Standard_Real anoffset = 0.;
3699 Standard_Real anAdjustPar = AdjustPeriodic(aParameter, alowerboundary, aupperboundary, aPeriod, anoffset);
3700
3701 Standard_Real adist = (bIsFirstBoundary) ? fabs(anAdjustPar - alowerboundary) : fabs(anAdjustPar - aupperboundary);
3702 Standard_Real anotherPar = (bIsFirstBoundary) ? (aupperboundary - adist) : (alowerboundary + adist);
3703 anotherPar += anoffset;
3704 Standard_Integer aneighbourpointindex = (j == 0) ? aListOfIndex.First() : aListOfIndex.Last();
3705 const IntSurf_PntOn2S& aNeighbourPoint = theWLine->Point(aneighbourpointindex);
3706 Standard_Real nU1, nV1;
3707
3708 if(surfit == 0)
3709 aNeighbourPoint.ParametersOnS1(nU1, nV1);
3710 else
3711 aNeighbourPoint.ParametersOnS2(nU1, nV1);
3712
3713 Standard_Real adist1 = (bIsUBoundary) ? fabs(nU1 - U) : fabs(nV1 - V);
3714 Standard_Real adist2 = (bIsUBoundary) ? fabs(nU1 - anotherPar) : fabs(nV1 - anotherPar);
3715 bComputeLineEnd = Standard_True;
3716 Standard_Boolean bCheckAngle1 = Standard_False;
3717 Standard_Boolean bCheckAngle2 = Standard_False;
3718 gp_Vec2d aNewVec;
3719 Standard_Real anewU = (bIsUBoundary) ? anotherPar : U;
3720 Standard_Real anewV = (bIsUBoundary) ? V : anotherPar;
3721
3722 if(((adist1 - adist2) > Precision::PConfusion()) &&
3723 (adist2 < (aPeriod / 4.))) {
3724 bCheckAngle1 = Standard_True;
3725 aNewVec = gp_Vec2d(gp_Pnt2d(nU1, nV1), gp_Pnt2d(anewU, anewV));
3726
3727 if(aNewVec.SquareMagnitude() < (gp::Resolution() * gp::Resolution())) {
3728 aNewP.SetValue((surfit == 0), anewU, anewV);
3729 bCheckAngle1 = Standard_False;
3730 }
3731 }
3732 else if(adist1 < (aPeriod / 4.)) {
3733 bCheckAngle2 = Standard_True;
3734 aNewVec = gp_Vec2d(gp_Pnt2d(nU1, nV1), gp_Pnt2d(U, V));
3735
3736 if(aNewVec.SquareMagnitude() < (gp::Resolution() * gp::Resolution())) {
3737 bCheckAngle2 = Standard_False;
3738 }
3739 }
3740
3741 if(bCheckAngle1 || bCheckAngle2) {
3742 // assume there are at least two points in line (see "if" above)
3743 Standard_Integer anindexother = aneighbourpointindex;
3744
3745 while((anindexother <= aListOfIndex.Last()) && (anindexother >= aListOfIndex.First())) {
3746 anindexother = (j == 0) ? (anindexother + 1) : (anindexother - 1);
3747 const IntSurf_PntOn2S& aPrevNeighbourPoint = theWLine->Point(anindexother);
3748 Standard_Real nU2, nV2;
3749
3750 if(surfit == 0)
3751 aPrevNeighbourPoint.ParametersOnS1(nU2, nV2);
3752 else
3753 aPrevNeighbourPoint.ParametersOnS2(nU2, nV2);
3754 gp_Vec2d aVecOld(gp_Pnt2d(nU2, nV2), gp_Pnt2d(nU1, nV1));
3755
3756 if(aVecOld.SquareMagnitude() <= (gp::Resolution() * gp::Resolution())) {
3757 continue;
3758 }
3759 else {
3760 Standard_Real anAngle = aNewVec.Angle(aVecOld);
3761
3762 if((fabs(anAngle) < (M_PI * 0.25)) && (aNewVec.Dot(aVecOld) > 0.)) {
3763
3764 if(bCheckAngle1) {
3765 Standard_Real U1, U2, V1, V2;
3766 IntSurf_PntOn2S atmppoint = aNewP;
3767 atmppoint.SetValue((surfit == 0), anewU, anewV);
3768 atmppoint.Parameters(U1, V1, U2, V2);
3769 gp_Pnt P1 = theSurface1->Value(U1, V1);
3770 gp_Pnt P2 = theSurface2->Value(U2, V2);
3771 gp_Pnt P0 = aPoint.Value();
3772
3773 if(P0.IsEqual(P1, aTol) &&
3774 P0.IsEqual(P2, aTol) &&
3775 P1.IsEqual(P2, aTol)) {
3776 bComputeLineEnd = Standard_False;
3777 aNewP.SetValue((surfit == 0), anewU, anewV);
3778 }
3779 }
3780
3781 if(bCheckAngle2) {
3782 bComputeLineEnd = Standard_False;
3783 }
3784 }
3785 break;
3786 }
3787 } // end while(anindexother...)
3788 }
3789 }
3790 }
3791 else if ( bIsNearBoundary ) {
3792 bComputeLineEnd = Standard_True;
3793 }
3794
3795 if(bComputeLineEnd) {
3796
3797 gp_Pnt2d anewpoint;
3798 Standard_Boolean found = Standard_False;
3799
3800 if ( bIsNearBoundary ) {
3801 // re-compute point near natural boundary or near tangent zone
3802 Standard_Real u1, v1, u2, v2;
3803 aNewP.Parameters( u1, v1, u2, v2 );
3804 if(surfit == 0)
3805 anewpoint = gp_Pnt2d( u1, v1 );
3806 else
3807 anewpoint = gp_Pnt2d( u2, v2 );
3808
3809 Standard_Integer aneighbourpointindex1 = (j == 0) ? aListOfIndex.First() : aListOfIndex.Last();
3810 const IntSurf_PntOn2S& aNeighbourPoint = theWLine->Point(aneighbourpointindex1);
3811 Standard_Real nU1, nV1;
3812
3813 if(surfit == 0)
3814 aNeighbourPoint.ParametersOnS1(nU1, nV1);
3815 else
3816 aNeighbourPoint.ParametersOnS2(nU1, nV1);
3817 gp_Pnt2d ap1(nU1, nV1);
3818 gp_Pnt2d ap2;
3819
3820
3821 if ( aZoneIndex ) {
3822 // exclude point from a tangent zone
3823 anewpoint = AdjustByNeighbour( ap1, anewpoint, aGASurface );
3824 gp_Pnt2d aPZone = (surfit == 0) ? aTanZoneS1->Value(aZoneIndex) : aTanZoneS2->Value(aZoneIndex);
3825 Standard_Real aZoneRadius = aTanZoneRadius->Value(aZoneIndex);
3826
3827 if ( FindPoint(ap1, anewpoint, umin, umax, vmin, vmax,
3828 aPZone, aZoneRadius, aGASurface, ap2) ) {
3829 anewpoint = ap2;
3830 found = Standard_True;
3831 }
3832 }
3833 else if ( aGASurface->IsUPeriodic() || aGASurface->IsVPeriodic() ) {
3834 // re-compute point near boundary if shifted on a period
3835 ap2 = AdjustByNeighbour( ap1, anewpoint, aGASurface );
3836
3837 if ( ( ap2.X() < umin ) || ( ap2.X() > umax ) ||
3838 ( ap2.Y() < vmin ) || ( ap2.Y() > vmax ) ) {
3839 found = FindPoint(ap1, ap2, umin, umax, vmin, vmax, anewpoint);
3840 }
3841 else {
3842 anewpoint = ap2;
3843 aNewP.SetValue( (surfit == 0), anewpoint.X(), anewpoint.Y() );
3844 }
3845 }
3846 }
3847 else {
3848
3849 Standard_Integer aneighbourpointindex1 = (j == 0) ? aListOfIndex.First() : aListOfIndex.Last();
3850 const IntSurf_PntOn2S& aNeighbourPoint = theWLine->Point(aneighbourpointindex1);
3851 Standard_Real nU1, nV1;
3852
3853 if(surfit == 0)
3854 aNeighbourPoint.ParametersOnS1(nU1, nV1);
3855 else
3856 aNeighbourPoint.ParametersOnS2(nU1, nV1);
3857 gp_Pnt2d ap1(nU1, nV1);
3858 gp_Pnt2d ap2(nU1, nV1);
3859 Standard_Integer aneighbourpointindex2 = aneighbourpointindex1;
3860
3861 while((aneighbourpointindex2 <= aListOfIndex.Last()) && (aneighbourpointindex2 >= aListOfIndex.First())) {
3862 aneighbourpointindex2 = (j == 0) ? (aneighbourpointindex2 + 1) : (aneighbourpointindex2 - 1);
3863 const IntSurf_PntOn2S& aPrevNeighbourPoint = theWLine->Point(aneighbourpointindex2);
3864 Standard_Real nU2, nV2;
3865
3866 if(surfit == 0)
3867 aPrevNeighbourPoint.ParametersOnS1(nU2, nV2);
3868 else
3869 aPrevNeighbourPoint.ParametersOnS2(nU2, nV2);
3870 ap2.SetX(nU2);
3871 ap2.SetY(nV2);
3872
3873 if(ap1.SquareDistance(ap2) > (gp::Resolution() * gp::Resolution())) {
3874 break;
3875 }
3876 }
3877 found = FindPoint(ap2, ap1, umin, umax, vmin, vmax, anewpoint);
3878 }
3879
3880 if(found) {
3881 // check point
3882 Standard_Real aCriteria = BRep_Tool::Tolerance(theFace1) + BRep_Tool::Tolerance(theFace2);
3883 GeomAPI_ProjectPointOnSurf& aProjector =
3884 (surfit == 0) ? aContext->ProjPS(theFace2) : aContext->ProjPS(theFace1);
3885 Handle(GeomAdaptor_HSurface) aSurface = (surfit == 0) ? theSurface1 : theSurface2;
3886
3887 Handle(GeomAdaptor_HSurface) aSurfaceOther = (surfit == 0) ? theSurface2 : theSurface1;
3888
3889 gp_Pnt aP3d = aSurface->Value(anewpoint.X(), anewpoint.Y());
3890 aProjector.Perform(aP3d);
3891
3892 if(aProjector.IsDone()) {
3893 if(aProjector.LowerDistance() < aCriteria) {
3894 Standard_Real foundU = U, foundV = V;
3895 aProjector.LowerDistanceParameters(foundU, foundV);
3896
3897 //Correction of projected coordinates. Begin
3898 //Note, it may be shifted on a period
3899 Standard_Integer aneindex1 = (j == 0) ? aListOfIndex.First() : aListOfIndex.Last();
3900 const IntSurf_PntOn2S& aNeighbourPoint = theWLine->Point(aneindex1);
3901 Standard_Real nUn, nVn;
3902
3903 if(surfit == 0)
3904 aNeighbourPoint.ParametersOnS2(nUn, nVn);
3905 else
3906 aNeighbourPoint.ParametersOnS1(nUn, nVn);
3907 gp_Pnt2d aNeighbour2d(nUn, nVn);
3908 gp_Pnt2d anAdjustedPoint = AdjustByNeighbour( aNeighbour2d, gp_Pnt2d(foundU, foundV), aSurfaceOther );
3909 foundU = anAdjustedPoint.X();
3910 foundV = anAdjustedPoint.Y();
3911
3912 if ( ( anAdjustedPoint.X() < umin ) && ( anAdjustedPoint.X() > umax ) &&
3913 ( anAdjustedPoint.Y() < vmin ) && ( anAdjustedPoint.Y() > vmax ) ) {
3914 // attempt to roughly re-compute point
3915 foundU = ( foundU < umin ) ? umin : foundU;
3916 foundU = ( foundU > umax ) ? umax : foundU;
3917 foundV = ( foundV < vmin ) ? vmin : foundV;
3918 foundV = ( foundV > vmax ) ? vmax : foundV;
3919
3920 GeomAPI_ProjectPointOnSurf& aProjector2 =
3921 (surfit == 0) ? aContext->ProjPS(theFace1) : aContext->ProjPS(theFace2);
3922
3923 aP3d = aSurfaceOther->Value(foundU, foundV);
3924 aProjector2.Perform(aP3d);
3925
3926 if(aProjector2.IsDone()) {
3927 if(aProjector2.LowerDistance() < aCriteria) {
3928 Standard_Real foundU2 = anewpoint.X(), foundV2 = anewpoint.Y();
3929 aProjector2.LowerDistanceParameters(foundU2, foundV2);
3930 anewpoint.SetX(foundU2);
3931 anewpoint.SetY(foundV2);
3932 }
3933 }
3934 }
3935 //Correction of projected coordinates. End
3936
3937 if(surfit == 0)
3938 aNewP.SetValue(aP3d, anewpoint.X(), anewpoint.Y(), foundU, foundV);
3939 else
3940 aNewP.SetValue(aP3d, foundU, foundV, anewpoint.X(), anewpoint.Y());
3941 }
3942 }
3943 }
3944 }
3945 }
3946 aSeqOfPntOn2S->Add(aNewP);
3947 aListOfFLIndex.Append(aSeqOfPntOn2S->NbPoints());
3948 }
3949 anArrayOfLineEnds.SetValue(i, aListOfFLIndex);
3950 }
3951 // Correct wlines.end
3952
3953 // Split wlines.begin
3954 Standard_Integer nbiter;
3955 //
3956 nbiter=1;
3957 if (!bAvoidLineConstructor) {
3958 nbiter=theLConstructor.NbParts();
3959 }
3960 //
3961 for(j = 1; j <= nbiter; ++j) {
3962 Standard_Real fprm, lprm;
3963 Standard_Integer ifprm, ilprm;
3964 //
3965 if(bAvoidLineConstructor) {
3966 ifprm = 1;
3967 ilprm = theWLine->NbPnts();
3968 }
3969 else {
3970 theLConstructor.Part(j, fprm, lprm);
3971 ifprm = (Standard_Integer)fprm;
3972 ilprm = (Standard_Integer)lprm;
3973 }
3974
3975 Handle(IntSurf_LineOn2S) aLineOn2S = new IntSurf_LineOn2S();
3976 //
3977 for(i = 1; i <= nblines; i++) {
3978 if(anArrayOfLineType.Value(i) != 0) {
3979 continue;
3980 }
3981 const TColStd_ListOfInteger& aListOfIndex = anArrayOfLines.Value(i);
3982
3983 if(aListOfIndex.Extent() < 2) {
3984 continue;
3985 }
3986 const TColStd_ListOfInteger& aListOfFLIndex = anArrayOfLineEnds.Value(i);
3987 Standard_Boolean bhasfirstpoint = (aListOfFLIndex.Extent() == 2);
3988 Standard_Boolean bhaslastpoint = (aListOfFLIndex.Extent() == 2);
3989
3990 if(!bhasfirstpoint && !aListOfFLIndex.IsEmpty()) {
3991 bhasfirstpoint = (i != 1);
3992 }
3993
3994 if(!bhaslastpoint && !aListOfFLIndex.IsEmpty()) {
3995 bhaslastpoint = (i != nblines);
3996 }
3997 Standard_Boolean bIsFirstInside = ((ifprm >= aListOfIndex.First()) && (ifprm <= aListOfIndex.Last()));
3998 Standard_Boolean bIsLastInside = ((ilprm >= aListOfIndex.First()) && (ilprm <= aListOfIndex.Last()));
3999
4000 if(!bIsFirstInside && !bIsLastInside) {
4001 if((ifprm < aListOfIndex.First()) && (ilprm > aListOfIndex.Last())) {
4002 // append whole line, and boundaries if neccesary
4003 if(bhasfirstpoint) {
4004 const IntSurf_PntOn2S& aP = aSeqOfPntOn2S->Value(aListOfFLIndex.First());
4005 aLineOn2S->Add(aP);
4006 }
4007 TColStd_ListIteratorOfListOfInteger anIt(aListOfIndex);
4008
4009 for(; anIt.More(); anIt.Next()) {
4010 const IntSurf_PntOn2S& aP = theWLine->Point(anIt.Value());
4011 aLineOn2S->Add(aP);
4012 }
4013
4014 if(bhaslastpoint) {
4015 const IntSurf_PntOn2S& aP = aSeqOfPntOn2S->Value(aListOfFLIndex.Last());
4016 aLineOn2S->Add(aP);
4017 }
4018
4019 // check end of split line (end is almost always)
4020 Standard_Integer aneighbour = i + 1;
4021 Standard_Boolean bIsEndOfLine = Standard_True;
4022
4023 if(aneighbour <= nblines) {
4024 const TColStd_ListOfInteger& aListOfNeighbourIndex = anArrayOfLines.Value(aneighbour);
4025
4026 if((anArrayOfLineType.Value(aneighbour) != 0) &&
4027 (aListOfNeighbourIndex.IsEmpty())) {
4028 bIsEndOfLine = Standard_False;
4029 }
4030 }
4031
4032 if(bIsEndOfLine) {
4033 if(aLineOn2S->NbPoints() > 1) {
4034 Handle(IntPatch_WLine) aNewWLine =
4035 new IntPatch_WLine(aLineOn2S, Standard_False);
4036 theNewLines.Append(aNewWLine);
4037 }
4038 aLineOn2S = new IntSurf_LineOn2S();
4039 }
4040 }
4041 continue;
4042 }
4043 // end if(!bIsFirstInside && !bIsLastInside)
4044
4045 if(bIsFirstInside && bIsLastInside) {
4046 // append inside points between ifprm and ilprm
4047 TColStd_ListIteratorOfListOfInteger anIt(aListOfIndex);
4048
4049 for(; anIt.More(); anIt.Next()) {
4050 if((anIt.Value() < ifprm) || (anIt.Value() > ilprm))
4051 continue;
4052 const IntSurf_PntOn2S& aP = theWLine->Point(anIt.Value());
4053 aLineOn2S->Add(aP);
4054 }
4055 }
4056 else {
4057
4058 if(bIsFirstInside) {
4059 // append points from ifprm to last point + boundary point
4060 TColStd_ListIteratorOfListOfInteger anIt(aListOfIndex);
4061
4062 for(; anIt.More(); anIt.Next()) {
4063 if(anIt.Value() < ifprm)
4064 continue;
4065 const IntSurf_PntOn2S& aP = theWLine->Point(anIt.Value());
4066 aLineOn2S->Add(aP);
4067 }
4068
4069 if(bhaslastpoint) {
4070 const IntSurf_PntOn2S& aP = aSeqOfPntOn2S->Value(aListOfFLIndex.Last());
4071 aLineOn2S->Add(aP);
4072 }
4073 // check end of split line (end is almost always)
4074 Standard_Integer aneighbour = i + 1;
4075 Standard_Boolean bIsEndOfLine = Standard_True;
4076
4077 if(aneighbour <= nblines) {
4078 const TColStd_ListOfInteger& aListOfNeighbourIndex = anArrayOfLines.Value(aneighbour);
4079
4080 if((anArrayOfLineType.Value(aneighbour) != 0) &&
4081 (aListOfNeighbourIndex.IsEmpty())) {
4082 bIsEndOfLine = Standard_False;
4083 }
4084 }
4085
4086 if(bIsEndOfLine) {
4087 if(aLineOn2S->NbPoints() > 1) {
4088 Handle(IntPatch_WLine) aNewWLine =
4089 new IntPatch_WLine(aLineOn2S, Standard_False);
4090 theNewLines.Append(aNewWLine);
4091 }
4092 aLineOn2S = new IntSurf_LineOn2S();
4093 }
4094 }
4095 // end if(bIsFirstInside)
4096
4097 if(bIsLastInside) {
4098 // append points from first boundary point to ilprm
4099 if(bhasfirstpoint) {
4100 const IntSurf_PntOn2S& aP = aSeqOfPntOn2S->Value(aListOfFLIndex.First());
4101 aLineOn2S->Add(aP);
4102 }
4103 TColStd_ListIteratorOfListOfInteger anIt(aListOfIndex);
4104
4105 for(; anIt.More(); anIt.Next()) {
4106 if(anIt.Value() > ilprm)
4107 continue;
4108 const IntSurf_PntOn2S& aP = theWLine->Point(anIt.Value());
4109 aLineOn2S->Add(aP);
4110 }
4111 }
4112 //end if(bIsLastInside)
4113 }
4114 }
4115
4116 if(aLineOn2S->NbPoints() > 1) {
4117 Handle(IntPatch_WLine) aNewWLine =
4118 new IntPatch_WLine(aLineOn2S, Standard_False);
4119 theNewLines.Append(aNewWLine);
4120 }
4121 }
4122 // Split wlines.end
4123
4124 return Standard_True;
4125}
4126
4127// ------------------------------------------------------------------------------------------------
4128// static function: ParameterOutOfBoundary
4129// purpose: Computes a new parameter for given curve. The corresponding 2d points
4130// does not lay on any boundary of given faces
4131// ------------------------------------------------------------------------------------------------
4132Standard_Boolean ParameterOutOfBoundary(const Standard_Real theParameter,
4133 const Handle(Geom_Curve)& theCurve,
4134 const TopoDS_Face& theFace1,
4135 const TopoDS_Face& theFace2,
4136 const Standard_Real theOtherParameter,
4137 const Standard_Boolean bIncreasePar,
4138 Standard_Real& theNewParameter,
4139 const Handle(IntTools_Context)& aContext)
4140{
4141 Standard_Boolean bIsComputed = Standard_False;
4142 theNewParameter = theParameter;
4143
4144 Standard_Real acurpar = theParameter;
4145 TopAbs_State aState = TopAbs_ON;
4146 Standard_Integer iter = 0;
4147 Standard_Real asumtol = BRep_Tool::Tolerance(theFace1) + BRep_Tool::Tolerance(theFace2);
4148 Standard_Real adelta = asumtol * 0.1;
4149 adelta = (adelta < Precision::Confusion()) ? Precision::Confusion() : adelta;
4150 Handle(Geom_Surface) aSurf1 = BRep_Tool::Surface(theFace1);
4151 Handle(Geom_Surface) aSurf2 = BRep_Tool::Surface(theFace2);
4152
4153 Standard_Real u1, u2, v1, v2;
4154
4155 GeomAPI_ProjectPointOnSurf aPrj1;
4156 aSurf1->Bounds(u1, u2, v1, v2);
4157 aPrj1.Init(aSurf1, u1, u2, v1, v2);
4158
4159 GeomAPI_ProjectPointOnSurf aPrj2;
4160 aSurf2->Bounds(u1, u2, v1, v2);
4161 aPrj2.Init(aSurf2, u1, u2, v1, v2);
4162
4163 while(aState == TopAbs_ON) {
4164 if(bIncreasePar)
4165 acurpar += adelta;
4166 else
4167 acurpar -= adelta;
4168 gp_Pnt aPCurrent = theCurve->Value(acurpar);
4169 aPrj1.Perform(aPCurrent);
4170 Standard_Real U=0., V=0.;
4171
4172 if(aPrj1.IsDone()) {
4173 aPrj1.LowerDistanceParameters(U, V);
4174 aState = aContext->StatePointFace(theFace1, gp_Pnt2d(U, V));
4175 }
4176
4177 if(aState != TopAbs_ON) {
4178 aPrj2.Perform(aPCurrent);
4179
4180 if(aPrj2.IsDone()) {
4181 aPrj2.LowerDistanceParameters(U, V);
4182 aState = aContext->StatePointFace(theFace2, gp_Pnt2d(U, V));
4183 }
4184 }
4185
4186 if(iter > 11) {
4187 break;
4188 }
4189 iter++;
4190 }
4191
4192 if(iter <= 11) {
4193 theNewParameter = acurpar;
4194 bIsComputed = Standard_True;
4195
4196 if(bIncreasePar) {
4197 if(acurpar >= theOtherParameter)
4198 theNewParameter = theOtherParameter;
4199 }
4200 else {
4201 if(acurpar <= theOtherParameter)
4202 theNewParameter = theOtherParameter;
4203 }
4204 }
4205 return bIsComputed;
4206}
4207
4208//=======================================================================
4209//function : IsCurveValid
4210//purpose :
4211//=======================================================================
4212Standard_Boolean IsCurveValid(Handle(Geom2d_Curve)& thePCurve)
4213{
4214 if(thePCurve.IsNull())
4215 return Standard_False;
4216
4217 Standard_Real tolint = 1.e-10;
4218 Geom2dAdaptor_Curve PCA;
4219 IntRes2d_Domain PCD;
4220 Geom2dInt_GInter PCI;
4221
4222 Standard_Real pf = 0., pl = 0.;
4223 gp_Pnt2d pntf, pntl;
4224
4225 if(!thePCurve->IsClosed() && !thePCurve->IsPeriodic()) {
4226 pf = thePCurve->FirstParameter();
4227 pl = thePCurve->LastParameter();
4228 pntf = thePCurve->Value(pf);
4229 pntl = thePCurve->Value(pl);
4230 PCA.Load(thePCurve);
4231 if(!PCA.IsPeriodic()) {
4232 if(PCA.FirstParameter() > pf) pf = PCA.FirstParameter();
4233 if(PCA.LastParameter() < pl) pl = PCA.LastParameter();
4234 }
4235 PCD.SetValues(pntf,pf,tolint,pntl,pl,tolint);
4236 PCI.Perform(PCA,PCD,tolint,tolint);
4237 if(PCI.IsDone())
4238 if(PCI.NbPoints() > 0) {
4239 return Standard_False;
4240 }
4241 }
4242
4243 return Standard_True;
4244}
4245
4246//=======================================================================
4247//static function : ApproxWithPCurves
4248//purpose : for bug 20964 only
4249//=======================================================================
4250Standard_Boolean ApproxWithPCurves(const gp_Cylinder& theCyl,
4251 const gp_Sphere& theSph)
4252{
4253 Standard_Boolean bRes = Standard_True;
4254 Standard_Real R1 = theCyl.Radius(), R2 = theSph.Radius();
4255
4256 if(R1 < 2.*R2) return bRes;
4257
4258 gp_Lin anCylAx(theCyl.Axis());
4259
4260 Standard_Real aDist = anCylAx.Distance(theSph.Location());
4261 Standard_Real aDRel = Abs(aDist - R1)/R2;
4262
4263 if(aDRel > .2) return bRes;
4264
4265 Standard_Real par = ElCLib::Parameter(anCylAx, theSph.Location());
4266 gp_Pnt aP = ElCLib::Value(par, anCylAx);
4267 gp_Vec aV(aP, theSph.Location());
4268
4269 Standard_Real dd = aV.Dot(theSph.Position().XDirection());
4270
4271 if(aDist < R1 && dd > 0.) return Standard_False;
4272 if(aDist > R1 && dd < 0.) return Standard_False;
4273
4274
4275 return bRes;
4276}
4277//=======================================================================
4278//function : PerformPlanes
4279//purpose :
4280//=======================================================================
4281void PerformPlanes(const Handle(GeomAdaptor_HSurface)& theS1,
4282 const Handle(GeomAdaptor_HSurface)& theS2,
4283 const Standard_Real TolAng,
4284 const Standard_Real TolTang,
4285 const Standard_Boolean theApprox1,
4286 const Standard_Boolean theApprox2,
4287 IntTools_SequenceOfCurves& theSeqOfCurve,
4288 Standard_Boolean& theTangentFaces)
4289{
4290
4291 gp_Pln aPln1 = theS1->Surface().Plane();
4292 gp_Pln aPln2 = theS2->Surface().Plane();
4293
4294 IntAna_QuadQuadGeo aPlnInter(aPln1, aPln2, TolAng, TolTang);
4295
4296 if(!aPlnInter.IsDone()) {
4297 theTangentFaces = Standard_False;
4298 return;
4299 }
4300
4301 IntAna_ResultType aResType = aPlnInter.TypeInter();
4302
4303 if(aResType == IntAna_Same) {
4304 theTangentFaces = Standard_True;
4305 return;
4306 }
4307
4308 theTangentFaces = Standard_False;
4309
4310 if(aResType == IntAna_Empty) {
4311 return;
4312 }
4313
4314 gp_Lin aLin = aPlnInter.Line(1);
4315
4316 ProjLib_Plane aProj;
4317
4318 aProj.Init(aPln1);
4319 aProj.Project(aLin);
4320 gp_Lin2d aLin2d1 = aProj.Line();
4321 //
4322 aProj.Init(aPln2);
4323 aProj.Project(aLin);
4324 gp_Lin2d aLin2d2 = aProj.Line();
4325 //
4326 //classify line2d1 relatively first plane
4327 Standard_Real P11, P12;
4328 Standard_Boolean IsCrossed = ClassifyLin2d(theS1, aLin2d1, TolTang, P11, P12);
4329 if(!IsCrossed) return;
4330 //classify line2d2 relatively second plane
4331 Standard_Real P21, P22;
4332 IsCrossed = ClassifyLin2d(theS2, aLin2d2, TolTang, P21, P22);
4333 if(!IsCrossed) return;
4334
4335 //Analysis of parametric intervals: must have common part
4336
4337 if(P21 >= P12) return;
4338 if(P22 <= P11) return;
4339
4340 Standard_Real pmin, pmax;
4341 pmin = Max(P11, P21);
4342 pmax = Min(P12, P22);
4343
4344 if(pmax - pmin <= TolTang) return;
4345
4346 Handle(Geom_Line) aGLin = new Geom_Line(aLin);
4347
4348 IntTools_Curve aCurve;
4349 Handle(Geom_TrimmedCurve) aGTLin = new Geom_TrimmedCurve(aGLin, pmin, pmax);
4350
4351 aCurve.SetCurve(aGTLin);
4352
4353 if(theApprox1) {
4354 Handle(Geom2d_Line) C2d = new Geom2d_Line(aLin2d1);
4355 aCurve.SetFirstCurve2d(new Geom2d_TrimmedCurve(C2d, pmin, pmax));
4356 }
4357 else {
4358 Handle(Geom2d_Curve) H1;
4359 aCurve.SetFirstCurve2d(H1);
4360 }
4361 if(theApprox2) {
4362 Handle(Geom2d_Line) C2d = new Geom2d_Line(aLin2d2);
4363 aCurve.SetSecondCurve2d(new Geom2d_TrimmedCurve(C2d, pmin, pmax));
4364 }
4365 else {
4366 Handle(Geom2d_Curve) H1;
4367 aCurve.SetFirstCurve2d(H1);
4368 }
4369
4370 theSeqOfCurve.Append(aCurve);
4371
4372}
4373
4374//=======================================================================
4375//function : ClassifyLin2d
4376//purpose :
4377//=======================================================================
4378static inline Standard_Boolean INTER(const Standard_Real d1,
4379 const Standard_Real d2,
4380 const Standard_Real tol)
4381{
4382 return (d1 > tol && d2 < -tol) ||
4383 (d1 < -tol && d2 > tol) ||
4384 ((d1 <= tol && d1 >= -tol) && (d2 > tol || d2 < -tol)) ||
4385 ((d2 <= tol && d2 >= -tol) && (d1 > tol || d1 < -tol));
4386}
4387static inline Standard_Boolean COINC(const Standard_Real d1,
4388 const Standard_Real d2,
4389 const Standard_Real tol)
4390{
4391 return (d1 <= tol && d1 >= -tol) && (d2 <= tol && d2 >= -tol);
4392}
4393Standard_Boolean ClassifyLin2d(const Handle(GeomAdaptor_HSurface)& theS,
4394 const gp_Lin2d& theLin2d,
4395 const Standard_Real theTol,
4396 Standard_Real& theP1,
4397 Standard_Real& theP2)
4398
4399{
4400 Standard_Real xmin, xmax, ymin, ymax, d1, d2, A, B, C;
4401 Standard_Real par[2];
4402 Standard_Integer nbi = 0;
4403
4404 xmin = theS->Surface().FirstUParameter();
4405 xmax = theS->Surface().LastUParameter();
4406 ymin = theS->Surface().FirstVParameter();
4407 ymax = theS->Surface().LastVParameter();
4408
4409 theLin2d.Coefficients(A, B, C);
4410
4411 //xmin, ymin <-> xmin, ymax
4412 d1 = A*xmin + B*ymin + C;
4413 d2 = A*xmin + B*ymax + C;
4414
4415 if(INTER(d1, d2, theTol)) {
4416 //Intersection with boundary
4417 Standard_Real y = -(C + A*xmin)/B;
4418 par[nbi] = ElCLib::Parameter(theLin2d, gp_Pnt2d(xmin, y));
4419 nbi++;
4420 }
4421 else if (COINC(d1, d2, theTol)) {
4422 //Coincidence with boundary
4423 par[0] = ElCLib::Parameter(theLin2d, gp_Pnt2d(xmin, ymin));
4424 par[1] = ElCLib::Parameter(theLin2d, gp_Pnt2d(xmin, ymax));
4425 nbi = 2;
4426 }
4427
4428 if(nbi == 2) {
4429
4430 if(fabs(par[0]-par[1]) > theTol) {
4431 theP1 = Min(par[0], par[1]);
4432 theP2 = Max(par[0], par[1]);
4433 return Standard_True;
4434 }
4435 else return Standard_False;
4436
4437 }
4438
4439 //xmin, ymax <-> xmax, ymax
4440 d1 = d2;
4441 d2 = A*xmax + B*ymax + C;
4442
4443 if(d1 > theTol || d1 < -theTol) {//to avoid checking of
4444 //coincidence with the same point
4445 if(INTER(d1, d2, theTol)) {
4446 Standard_Real x = -(C + B*ymax)/A;
4447 par[nbi] = ElCLib::Parameter(theLin2d, gp_Pnt2d(x, ymax));
4448 nbi++;
4449 }
4450 else if (COINC(d1, d2, theTol)) {
4451 par[0] = ElCLib::Parameter(theLin2d, gp_Pnt2d(xmin, ymax));
4452 par[1] = ElCLib::Parameter(theLin2d, gp_Pnt2d(xmax, ymax));
4453 nbi = 2;
4454 }
4455 }
4456
4457 if(nbi == 2) {
4458
4459 if(fabs(par[0]-par[1]) > theTol) {
4460 theP1 = Min(par[0], par[1]);
4461 theP2 = Max(par[0], par[1]);
4462 return Standard_True;
4463 }
4464 else return Standard_False;
4465
4466 }
4467
4468 //xmax, ymax <-> xmax, ymin
4469 d1 = d2;
4470 d2 = A*xmax + B*ymin + C;
4471
4472 if(d1 > theTol || d1 < -theTol) {
4473 if(INTER(d1, d2, theTol)) {
4474 Standard_Real y = -(C + A*xmax)/B;
4475 par[nbi] = ElCLib::Parameter(theLin2d, gp_Pnt2d(xmax, y));
4476 nbi++;
4477 }
4478 else if (COINC(d1, d2, theTol)) {
4479 par[0] = ElCLib::Parameter(theLin2d, gp_Pnt2d(xmax, ymax));
4480 par[1] = ElCLib::Parameter(theLin2d, gp_Pnt2d(xmax, ymin));
4481 nbi = 2;
4482 }
4483 }
4484
4485 if(nbi == 2) {
4486 if(fabs(par[0]-par[1]) > theTol) {
4487 theP1 = Min(par[0], par[1]);
4488 theP2 = Max(par[0], par[1]);
4489 return Standard_True;
4490 }
4491 else return Standard_False;
4492 }
4493
4494 //xmax, ymin <-> xmin, ymin
4495 d1 = d2;
4496 d2 = A*xmin + B*ymin + C;
4497
4498 if(d1 > theTol || d1 < -theTol) {
4499 if(INTER(d1, d2, theTol)) {
4500 Standard_Real x = -(C + B*ymin)/A;
4501 par[nbi] = ElCLib::Parameter(theLin2d, gp_Pnt2d(x, ymin));
4502 nbi++;
4503 }
4504 else if (COINC(d1, d2, theTol)) {
4505 par[0] = ElCLib::Parameter(theLin2d, gp_Pnt2d(xmax, ymin));
4506 par[1] = ElCLib::Parameter(theLin2d, gp_Pnt2d(xmin, ymin));
4507 nbi = 2;
4508 }
4509 }
4510
4511 if(nbi == 2) {
4512 if(fabs(par[0]-par[1]) > theTol) {
4513 theP1 = Min(par[0], par[1]);
4514 theP2 = Max(par[0], par[1]);
4515 return Standard_True;
4516 }
4517 else return Standard_False;
4518 }
4519
4520 return Standard_False;
4521
4522}
4523//
4524//=======================================================================
4525//function : ApproxParameters
4526//purpose :
4527//=======================================================================
4528void ApproxParameters(const Handle(GeomAdaptor_HSurface)& aHS1,
4529 const Handle(GeomAdaptor_HSurface)& aHS2,
4530 Standard_Integer& iDegMin,
4531 Standard_Integer& iDegMax,
4532 Standard_Integer& iNbIter)
4533
4534{
4535 GeomAbs_SurfaceType aTS1, aTS2;
4536
4537 //
4538 iNbIter=0;
4539 iDegMin=4;
4540 iDegMax=8;
4541 //
4542 aTS1=aHS1->Surface().GetType();
4543 aTS2=aHS2->Surface().GetType();
4544 //
4545 // Cylinder/Torus
4546 if ((aTS1==GeomAbs_Cylinder && aTS2==GeomAbs_Torus) ||
4547 (aTS2==GeomAbs_Cylinder && aTS1==GeomAbs_Torus)) {
4548 Standard_Real aRC, aRT, dR, aPC;
4549 gp_Cylinder aCylinder;
4550 gp_Torus aTorus;
4551 //
4552 aPC=Precision::Confusion();
4553 //
4554 aCylinder=(aTS1==GeomAbs_Cylinder)? aHS1->Surface().Cylinder() : aHS2->Surface().Cylinder();
4555 aTorus=(aTS1==GeomAbs_Torus)? aHS1->Surface().Torus() : aHS2->Surface().Torus();
4556 //
4557 aRC=aCylinder.Radius();
4558 aRT=aTorus.MinorRadius();
4559 dR=aRC-aRT;
4560 if (dR<0.) {
4561 dR=-dR;
4562 }
4563 //
4564 if (dR<aPC) {
4565 iDegMax=6;
4566 }
4567 }
4568 if (aTS1==GeomAbs_Cylinder && aTS2==GeomAbs_Cylinder) {
4569 iNbIter=1;
4570 }
4571}
4572//=======================================================================
4573//function : Tolerances
4574//purpose :
4575//=======================================================================
4576void Tolerances(const Handle(GeomAdaptor_HSurface)& aHS1,
4577 const Handle(GeomAdaptor_HSurface)& aHS2,
4578 Standard_Real& ,//aTolArc,
4579 Standard_Real& aTolTang,
4580 Standard_Real& ,//aUVMaxStep,
4581 Standard_Real& )//aDeflection)
4582{
4583 GeomAbs_SurfaceType aTS1, aTS2;
4584 //
4585 aTS1=aHS1->Surface().GetType();
4586 aTS2=aHS2->Surface().GetType();
4587 //
4588 // Cylinder/Torus
4589 if ((aTS1==GeomAbs_Cylinder && aTS2==GeomAbs_Torus) ||
4590 (aTS2==GeomAbs_Cylinder && aTS1==GeomAbs_Torus)) {
4591 Standard_Real aRC, aRT, dR, aPC;
4592 gp_Cylinder aCylinder;
4593 gp_Torus aTorus;
4594 //
4595 aPC=Precision::Confusion();
4596 //
4597 aCylinder=(aTS1==GeomAbs_Cylinder)? aHS1->Surface().Cylinder() : aHS2->Surface().Cylinder();
4598 aTorus=(aTS1==GeomAbs_Torus)? aHS1->Surface().Torus() : aHS2->Surface().Torus();
4599 //
4600 aRC=aCylinder.Radius();
4601 aRT=aTorus.MinorRadius();
4602 dR=aRC-aRT;
4603 if (dR<0.) {
4604 dR=-dR;
4605 }
4606 //
4607 if (dR<aPC) {
4608 aTolTang=0.1*aTolTang;
4609 }
4610 }
4611}
4612//=======================================================================
4613//function : SortTypes
4614//purpose :
4615//=======================================================================
4616Standard_Boolean SortTypes(const GeomAbs_SurfaceType aType1,
4617 const GeomAbs_SurfaceType aType2)
4618{
4619 Standard_Boolean bRet;
4620 Standard_Integer aI1, aI2;
4621 //
4622 bRet=Standard_False;
4623 //
4624 aI1=IndexType(aType1);
4625 aI2=IndexType(aType2);
4626 if (aI1<aI2){
4627 bRet=!bRet;
4628 }
4629 return bRet;
4630}
4631//=======================================================================
4632//function : IndexType
4633//purpose :
4634//=======================================================================
4635Standard_Integer IndexType(const GeomAbs_SurfaceType aType)
4636{
4637 Standard_Integer aIndex;
4638 //
4639 aIndex=11;
4640 //
4641 if (aType==GeomAbs_Plane) {
4642 aIndex=0;
4643 }
4644 else if (aType==GeomAbs_Cylinder) {
4645 aIndex=1;
4646 }
4647 else if (aType==GeomAbs_Cone) {
4648 aIndex=2;
4649 }
4650 else if (aType==GeomAbs_Sphere) {
4651 aIndex=3;
4652 }
4653 else if (aType==GeomAbs_Torus) {
4654 aIndex=4;
4655 }
4656 else if (aType==GeomAbs_BezierSurface) {
4657 aIndex=5;
4658 }
4659 else if (aType==GeomAbs_BSplineSurface) {
4660 aIndex=6;
4661 }
4662 else if (aType==GeomAbs_SurfaceOfRevolution) {
4663 aIndex=7;
4664 }
4665 else if (aType==GeomAbs_SurfaceOfExtrusion) {
4666 aIndex=8;
4667 }
4668 else if (aType==GeomAbs_OffsetSurface) {
4669 aIndex=9;
4670 }
4671 else if (aType==GeomAbs_OtherSurface) {
4672 aIndex=10;
4673 }
4674 return aIndex;
4675}
4676//=======================================================================
4677//function : DumpWLine
4678//purpose :
4679//=======================================================================
4680void DumpWLine(const Handle(IntPatch_WLine)& aWLine)
4681{
4682 Standard_Integer i, aNbPnts;
4683 Standard_Real aX, aY, aZ, aU1, aV1, aU2, aV2;
4684 //
4685 printf(" *WLine\n");
4686 aNbPnts=aWLine->NbPnts();
4687 for (i=1; i<=aNbPnts; ++i) {
4688 const IntSurf_PntOn2S aPntOn2S=aWLine->Point(i);
4689 const gp_Pnt& aP3D=aPntOn2S.Value();
4690 aP3D.Coord(aX, aY, aZ);
4691 aPntOn2S.Parameters(aU1, aV1, aU2, aV2);
4692 //
4693 printf("point p_%d %lf %lf %lf\n", i, aX, aY, aZ);
4694 //printf("point p_%d %20.15lf %20.15lf %20.15lf %20.15lf %20.15lf %20.15lf %20.15lf\n",
4695 // i, aX, aY, aZ, aU1, aV1, aU2, aV2);
4696 }
4697}
4698//=======================================================================
4699//function : RefineVector
4700//purpose :
4701//=======================================================================
4702void RefineVector(gp_Vec2d& aV2D)
4703{
4704 Standard_Integer k,m;
4705 Standard_Real aC[2], aEps, aR1, aR2, aNum;
4706 //
4707 aEps=RealEpsilon();
4708 aR1=1.-aEps;
4709 aR2=1.+aEps;
4710 //
4711 aV2D.Coord(aC[0], aC[1]);
4712 //
4713 for (k=0; k<2; ++k) {
4714 m=(k+1)%2;
4715 aNum=fabs(aC[k]);
4716 if (aNum>aR1 && aNum<aR2) {
4717 if (aC[k]<0.) {
4718 aC[k]=-1.;
4719 }
4720 else {
4721 aC[k]=1.;
4722 }
4723 aC[m]=0.;
4724 break;
4725 }
4726 }
4727 aV2D.SetCoord(aC[0], aC[1]);
4728}
4729//=======================================================================
4730//function : FindMaxSquareDistance
4731//purpose :
4732//=======================================================================
4733Standard_Real FindMaxSquareDistance (const Standard_Real aT1,
4734 const Standard_Real aT2,
4735 const Standard_Real aEps,
4736 const Handle(Geom_Curve)& aC3D,
4737 const Handle(Geom2d_Curve)& aC2D1,
4738 const Handle(Geom2d_Curve)& aC2D2,
4739 const Handle(GeomAdaptor_HSurface)& myHS1,
4740 const Handle(GeomAdaptor_HSurface)& myHS2,
4741 const TopoDS_Face& myFace1,
4742 const TopoDS_Face& myFace2,
4743 const Handle(IntTools_Context)& myContext)
4744{
4745 Standard_Real aA, aB, aCf, aX1, aX2, aF1, aF2, aX, aF;
4746 //
4747 aCf=1.6180339887498948482045868343656;// =0.5*(1.+sqrt(5.));
4748 aA=aT1;
4749 aB=aT2;
4750 aX1=aB-(aB-aA)/aCf;
4751 aF1=MaxSquareDistance(aX1,
4752 aC3D, aC2D1, aC2D2, myHS1, myHS2, myFace1, myFace2, myContext);
4753 aX2=aA+(aB-aA)/aCf;
4754 aF2=MaxSquareDistance(aX2,
4755 aC3D, aC2D1, aC2D2, myHS1, myHS2, myFace1, myFace2, myContext);
4756 //
4757 while(1) {
4758 //
4759 if (fabs(aA-aB)<aEps) {
4760 aX=0.5*(aA+aB);
4761 aF=MaxSquareDistance(aX,
4762 aC3D, aC2D1, aC2D2, myHS1, myHS2, myFace1, myFace2, myContext);
4763 break;
4764 }
4765 if (aF1<aF2){
4766 aA=aX1;
4767 aX1=aX2;
4768 aF1=aF2;
4769 aX2=aA+(aB-aA)/aCf;
4770 aF2=MaxSquareDistance(aX2,
4771 aC3D, aC2D1, aC2D2, myHS1, myHS2, myFace1, myFace2, myContext);
4772
4773 }
4774 else {
4775 aB=aX2;
4776 aX2=aX1;
4777 aF2=aF1;
4778 aX1=aB-(aB-aA)/aCf;
4779 aF1=MaxSquareDistance(aX1,
4780 aC3D, aC2D1, aC2D2, myHS1, myHS2, myFace1, myFace2, myContext);
4781 }
4782 }
4783 return aF;
4784}
4785//=======================================================================
4786//function : MaxSquareDistance
4787//purpose :
4788//=======================================================================
4789Standard_Real MaxSquareDistance (const Standard_Real aT,
4790 const Handle(Geom_Curve)& aC3D,
4791 const Handle(Geom2d_Curve)& aC2D1,
4792 const Handle(Geom2d_Curve)& aC2D2,
4793 const Handle(GeomAdaptor_HSurface) myHS1,
4794 const Handle(GeomAdaptor_HSurface) myHS2,
4795 const TopoDS_Face& aF1,
4796 const TopoDS_Face& aF2,
4797 const Handle(IntTools_Context)& aCtx)
4798{
4799 Standard_Boolean bIsDone;
4800 Standard_Integer i;
4801 Standard_Real aU, aV, aD2Max, aD2;
4802 gp_Pnt2d aP2D;
4803 gp_Pnt aP, aPS;
4804 //
4805 aD2Max=0.;
4806 //
4807 aC3D->D0(aT, aP);
4808 if (aC3D.IsNull()) {
4809 return aD2Max;
4810 }
4811 //
4812 for (i=0; i<2; ++i) {
4813 const Handle(GeomAdaptor_HSurface)& aGHS=(!i) ? myHS1 : myHS2;
4814 const TopoDS_Face &aF=(!i) ? aF1 : aF2;
4815 const Handle(Geom2d_Curve)& aC2D=(!i) ? aC2D1 : aC2D2;
4816 //
4817 if (!aC2D.IsNull()) {
4818 aC2D->D0(aT, aP2D);
4819 aP2D.Coord(aU, aV);
4820 aGHS->D0(aU, aV, aPS);
4821 aD2=aP.SquareDistance(aPS);
4822 if (aD2>aD2Max) {
4823 aD2Max=aD2;
4824 }
4825 }
4826 //
4827 GeomAPI_ProjectPointOnSurf& aProjector=aCtx->ProjPS(aF);
4828 //
4829 aProjector.Perform(aP);
4830 bIsDone=aProjector.IsDone();
4831 if (bIsDone) {
4832 aProjector.LowerDistanceParameters(aU, aV);
4833 aGHS->D0(aU, aV, aPS);
4834 aD2=aP.SquareDistance(aPS);
4835 if (aD2>aD2Max) {
4836 aD2Max=aD2;
4837 }
4838 }
4839 }
4840 //
4841 return aD2Max;
4842}