1 // Created on: 1999-03-05
2 // Created by: Fabrice SERVANT
3 // Copyright (c) 1999-1999 Matra Datavision
4 // Copyright (c) 1999-2014 OPEN CASCADE SAS
5 //
6 // This file is part of Open CASCADE Technology software library.
7 //
8 // This library is free software; you can redistribute it and/or modify it under
9 // the terms of the GNU Lesser General Public License version 2.1 as published
10 // by the Free Software Foundation, with special exception defined in the file
11 // OCCT_LGPL_EXCEPTION.txt. Consult the file LICENSE_LGPL_21.txt included in OCCT
12 // distribution for complete text of the license and disclaimer of any warranty.
13 //
14 // Alternatively, this file may be used under the terms of Open CASCADE
15 // commercial license or contractual agreement.
16
17 // modified by Edward AGAPOV (eap) Tue Jan 22 2002 (bug occ53)
18 // - improve SectionLine table management (avoid memory reallocation)
19 // - some protection against arrays overflow
20 // modified by Edward AGAPOV (eap) Thu Feb 14 2002 (occ139)
21 // - make Section Line parts rightly connected (prepend 2nd part to the 1st)
22 // - TriangleShape() for debugging purpose
23 // Modified by skv - Thu Sep 25 17:42:42 2003 OCC567
24 // modified by ofv Thu Apr 8 14:58:13 2004 fip
25
26 #include <Adaptor3d_Surface.hxx>
27 #include <Bnd_BoundSortBox.hxx>
28 #include <Bnd_Box.hxx>
29 #include <Bnd_HArray1OfBox.hxx>
30 #include <Bnd_Tools.hxx>
31 #include <BVH_BoxSet.hxx>
32 #include <BVH_LinearBuilder.hxx>
33 #include <BVH_Traverse.hxx>
34 #include <gp.hxx>
35 #include <gp_Pnt.hxx>
36 #include <IntPolyh_ListOfCouples.hxx>
37 #include <IntPolyh_Couple.hxx>
38 #include <IntPolyh_Edge.hxx>
39 #include <IntPolyh_MaillageAffinage.hxx>
40 #include <IntPolyh_Point.hxx>
41 #include <IntPolyh_SectionLine.hxx>
42 #include <IntPolyh_StartPoint.hxx>
43 #include <IntPolyh_Tools.hxx>
44 #include <IntPolyh_Triangle.hxx>
45 #include <Precision.hxx>
46 #include <TColStd_Array1OfInteger.hxx>
47 #include <TColStd_MapOfInteger.hxx>
48 #include <TColStd_ListIteratorOfListOfInteger.hxx>
49 #include <algorithm>
50 #include <NCollection_IndexedDataMap.hxx>
51
52 typedef NCollection_Array1<Standard_Integer> IntPolyh_ArrayOfInteger;
53 typedef NCollection_IndexedDataMap
54 <Standard_Integer,
55 TColStd_ListOfInteger,
56 TColStd_MapIntegerHasher> IntPolyh_IndexedDataMapOfIntegerListOfInteger;
57
58
59 static Standard_Real MyTolerance=10.0e-7;
60 static Standard_Real MyConfusionPrecision=10.0e-12;
61 static Standard_Real SquareMyConfusionPrecision=10.0e-24;
62 //
63 static
64 inline Standard_Real maxSR(const Standard_Real a,
65 const Standard_Real b,
66 const Standard_Real c);
67
68 static
69 inline Standard_Real minSR(const Standard_Real a,
70 const Standard_Real b,
71 const Standard_Real c);
72 static
73 Standard_Integer project6(const IntPolyh_Point &ax,
74 const IntPolyh_Point &p1,
75 const IntPolyh_Point &p2,
76 const IntPolyh_Point &p3,
77 const IntPolyh_Point &q1,
78 const IntPolyh_Point &q2,
79 const IntPolyh_Point &q3);
80 static
81 void TestNbPoints(const Standard_Integer ,
82 Standard_Integer &NbPoints,
83 Standard_Integer &NbPointsTotal,
84 const IntPolyh_StartPoint &Pt1,
85 const IntPolyh_StartPoint &Pt2,
86 IntPolyh_StartPoint &SP1,
87 IntPolyh_StartPoint &SP2);
88 static
89 void CalculPtsInterTriEdgeCoplanaires(const Standard_Integer TriSurfID,
90 const IntPolyh_Point &NormaleTri,
91 const IntPolyh_Triangle &Tri1,
92 const IntPolyh_Triangle &Tri2,
93 const IntPolyh_Point &PE1,
94 const IntPolyh_Point &PE2,
95 const IntPolyh_Point &Edge,
96 const Standard_Integer EdgeIndex,
97 const IntPolyh_Point &PT1,
98 const IntPolyh_Point &PT2,
99 const IntPolyh_Point &Cote,
100 const Standard_Integer CoteIndex,
101 IntPolyh_StartPoint &SP1,
102 IntPolyh_StartPoint &SP2,
103 Standard_Integer &NbPoints);
104 static
105 Standard_Boolean CheckCoupleAndGetAngle(const Standard_Integer T1,
106 const Standard_Integer T2,
107 Standard_Real& Angle,
108 IntPolyh_ListOfCouples &TTrianglesContacts);
109 static
110 Standard_Boolean CheckCoupleAndGetAngle2(const Standard_Integer T1,
111 const Standard_Integer T2,
112 const Standard_Integer T11,
113 const Standard_Integer T22,
114 IntPolyh_ListIteratorOfListOfCouples& theItCT11,
115 IntPolyh_ListIteratorOfListOfCouples& theItCT22,
116 Standard_Real & Angle,
117 IntPolyh_ListOfCouples &TTrianglesContacts);
118 static
119 Standard_Integer CheckNextStartPoint(IntPolyh_SectionLine & SectionLine,
120 IntPolyh_ArrayOfTangentZones & TTangentZones,
121 IntPolyh_StartPoint & SP,
122 const Standard_Boolean Prepend=Standard_False);
123
124 static
125 Standard_Boolean IsDegenerated(const Handle(Adaptor3d_Surface)& aS,
126 const Standard_Integer aIndex,
127 const Standard_Real aTol2,
128 Standard_Real& aDegX);
129 static
130 void DegeneratedIndex(const TColStd_Array1OfReal& Xpars,
131 const Standard_Integer aNbX,
132 const Handle(Adaptor3d_Surface)& aS,
133 const Standard_Integer aIsoDirection,
134 Standard_Integer& aI1,
135 Standard_Integer& aI2);
136
137 //=======================================================================
138 //class : IntPolyh_BoxBndTree
139 //purpose : BVH structure to contain the boxes of triangles
140 //=======================================================================
141 typedef BVH_BoxSet <Standard_Real, 3, Standard_Integer> IntPolyh_BoxBndTree;
142
143 //=======================================================================
144 //class : IntPolyh_BoxBndTreeSelector
145 //purpose : Selector of interfering boxes
146 //=======================================================================
147 class IntPolyh_BoxBndTreeSelector :
148 public BVH_PairTraverse<Standard_Real, 3, IntPolyh_BoxBndTree>
149 {
150 public:
151 typedef BVH_Box<Standard_Real, 3>::BVH_VecNt BVH_Vec3d;
152
153 //! Auxiliary structure to keep the pair of indices
154 struct PairIDs
155 {
PairIDsIntPolyh_BoxBndTreeSelector::PairIDs156 PairIDs (const Standard_Integer theId1 = -1,
157 const Standard_Integer theId2 = -1)
158 : ID1 (theId1), ID2 (theId2)
159 {}
160
operator <IntPolyh_BoxBndTreeSelector::PairIDs161 Standard_Boolean operator< (const PairIDs& theOther) const
162 {
163 return ID1 < theOther.ID1 ||
164 (ID1 == theOther.ID1 && ID2 < theOther.ID2);
165 }
166
167 Standard_Integer ID1;
168 Standard_Integer ID2;
169 };
170
171 public:
172
173 //! Constructor
IntPolyh_BoxBndTreeSelector()174 IntPolyh_BoxBndTreeSelector ()
175 {}
176
177 //! Rejects the node
RejectNode(const BVH_Vec3d & theCMin1,const BVH_Vec3d & theCMax1,const BVH_Vec3d & theCMin2,const BVH_Vec3d & theCMax2,Standard_Real &) const178 virtual Standard_Boolean RejectNode (const BVH_Vec3d& theCMin1,
179 const BVH_Vec3d& theCMax1,
180 const BVH_Vec3d& theCMin2,
181 const BVH_Vec3d& theCMax2,
182 Standard_Real&) const Standard_OVERRIDE
183 {
184 return BVH_Box<Standard_Real, 3> (theCMin1, theCMax1).IsOut (theCMin2, theCMax2);
185 }
186
187 //! Accepts the element
Accept(const Standard_Integer theID1,const Standard_Integer theID2)188 virtual Standard_Boolean Accept (const Standard_Integer theID1,
189 const Standard_Integer theID2) Standard_OVERRIDE
190 {
191 if (!myBVHSet1->Box (theID1).IsOut (myBVHSet2->Box (theID2)))
192 {
193 myPairs.push_back (PairIDs (myBVHSet1->Element (theID1), myBVHSet2->Element (theID2)));
194 return Standard_True;
195 }
196 return Standard_False;
197 }
198
199 //! Returns indices
Pairs() const200 const std::vector<PairIDs>& Pairs() const
201 {
202 return myPairs;
203 }
204
205 //! Sorts the resulting indices
Sort()206 void Sort()
207 {
208 std::sort (myPairs.begin(), myPairs.end());
209 }
210
211 private:
212
213 std::vector<PairIDs> myPairs;
214 };
215
216 //=======================================================================
217 //function : GetInterferingTriangles
218 //purpose : Returns indices of the triangles with interfering bounding boxes
219 //=======================================================================
220 static
GetInterferingTriangles(IntPolyh_ArrayOfTriangles & theTriangles1,const IntPolyh_ArrayOfPoints & thePoints1,IntPolyh_ArrayOfTriangles & theTriangles2,const IntPolyh_ArrayOfPoints & thePoints2,IntPolyh_IndexedDataMapOfIntegerListOfInteger & theCouples)221 void GetInterferingTriangles(IntPolyh_ArrayOfTriangles& theTriangles1,
222 const IntPolyh_ArrayOfPoints& thePoints1,
223 IntPolyh_ArrayOfTriangles& theTriangles2,
224 const IntPolyh_ArrayOfPoints& thePoints2,
225 IntPolyh_IndexedDataMapOfIntegerListOfInteger& theCouples)
226 {
227 // Use linear builder for BVH construction
228 opencascade::handle<BVH_LinearBuilder<Standard_Real, 3>> aLBuilder =
229 new BVH_LinearBuilder<Standard_Real, 3> (10);
230
231 // To find the triangles with interfering bounding boxes
232 // use the BVH structure
233 IntPolyh_BoxBndTree aBBTree1 (aLBuilder), aBBTree2 (aLBuilder);
234
235 // 1. Fill the trees with the boxes of the surfaces triangles
236 for (Standard_Integer i = 0; i < 2; ++i)
237 {
238 IntPolyh_BoxBndTree &aBBTree = !i ? aBBTree1 : aBBTree2;
239 IntPolyh_ArrayOfTriangles& aTriangles = !i ? theTriangles1 : theTriangles2;
240 const IntPolyh_ArrayOfPoints& aPoints = !i ? thePoints1 : thePoints2;
241
242 const Standard_Integer aNbT = aTriangles.NbItems();
243 aBBTree.SetSize (aNbT);
244 for (Standard_Integer j = 0; j < aNbT; ++j)
245 {
246 IntPolyh_Triangle& aT = aTriangles[j];
247 if (!aT.IsIntersectionPossible() || aT.IsDegenerated())
248 continue;
249
250 aBBTree.Add (j, Bnd_Tools::Bnd2BVH (aT.BoundingBox(aPoints)));
251 }
252
253 if (!aBBTree.Size())
254 return;
255 }
256 // 2. Construct BVH trees
257 aBBTree1.Build();
258 aBBTree2.Build();
259
260 // 3. Perform selection of the interfering triangles
261 IntPolyh_BoxBndTreeSelector aSelector;
262 aSelector.SetBVHSets (&aBBTree1, &aBBTree2);
263 aSelector.Select();
264 aSelector.Sort();
265
266 const std::vector<IntPolyh_BoxBndTreeSelector::PairIDs>& aPairs = aSelector.Pairs();
267 const Standard_Integer aNbPairs = static_cast<Standard_Integer>(aPairs.size());
268
269 for (Standard_Integer i = 0; i < aNbPairs; ++i)
270 {
271 const IntPolyh_BoxBndTreeSelector::PairIDs& aPair = aPairs[i];
272 TColStd_ListOfInteger* pTriangles2 = theCouples.ChangeSeek (aPair.ID1);
273 if (!pTriangles2)
274 pTriangles2 = &theCouples( theCouples.Add (aPair.ID1, TColStd_ListOfInteger()));
275 pTriangles2->Append (aPair.ID2);
276 }
277 }
278
279 //=======================================================================
280 //function : IntPolyh_MaillageAffinage
281 //purpose :
282 //=======================================================================
IntPolyh_MaillageAffinage(const Handle (Adaptor3d_Surface)& Surface1,const Handle (Adaptor3d_Surface)& Surface2,const Standard_Integer)283 IntPolyh_MaillageAffinage::IntPolyh_MaillageAffinage
284 (const Handle(Adaptor3d_Surface)& Surface1,
285 const Handle(Adaptor3d_Surface)& Surface2,
286 const Standard_Integer )
287 :
288 MaSurface1(Surface1),
289 MaSurface2(Surface2),
290 NbSamplesU1(10),
291 NbSamplesU2(10),
292 NbSamplesV1(10),
293 NbSamplesV2(10),
294 FlecheMax1(0.0),
295 FlecheMax2(0.0),
296 FlecheMin1(0.0),
297 FlecheMin2(0.0),
298 myEnlargeZone(Standard_False)
299 {
300 }
301 //=======================================================================
302 //function : IntPolyh_MaillageAffinage
303 //purpose :
304 //=======================================================================
IntPolyh_MaillageAffinage(const Handle (Adaptor3d_Surface)& Surface1,const Standard_Integer NbSU1,const Standard_Integer NbSV1,const Handle (Adaptor3d_Surface)& Surface2,const Standard_Integer NbSU2,const Standard_Integer NbSV2,const Standard_Integer)305 IntPolyh_MaillageAffinage::IntPolyh_MaillageAffinage
306 (const Handle(Adaptor3d_Surface)& Surface1,
307 const Standard_Integer NbSU1,
308 const Standard_Integer NbSV1,
309 const Handle(Adaptor3d_Surface)& Surface2,
310 const Standard_Integer NbSU2,
311 const Standard_Integer NbSV2,
312 const Standard_Integer )
313 :
314 MaSurface1(Surface1),
315 MaSurface2(Surface2),
316 NbSamplesU1(NbSU1),
317 NbSamplesU2(NbSU2),
318 NbSamplesV1(NbSV1),
319 NbSamplesV2(NbSV2),
320 FlecheMax1(0.0),
321 FlecheMax2(0.0),
322 FlecheMin1(0.0),
323 FlecheMin2(0.0),
324 myEnlargeZone(Standard_False)
325 {
326 }
327 //=======================================================================
328 //function : MakeSampling
329 //purpose :
330 //=======================================================================
MakeSampling(const Standard_Integer SurfID,TColStd_Array1OfReal & theUPars,TColStd_Array1OfReal & theVPars)331 void IntPolyh_MaillageAffinage::MakeSampling(const Standard_Integer SurfID,
332 TColStd_Array1OfReal& theUPars,
333 TColStd_Array1OfReal& theVPars)
334 {
335 if (SurfID == 1)
336 IntPolyh_Tools::MakeSampling(MaSurface1, NbSamplesU1, NbSamplesV1, myEnlargeZone, theUPars, theVPars);
337 else
338 IntPolyh_Tools::MakeSampling(MaSurface2, NbSamplesU2, NbSamplesV2, myEnlargeZone, theUPars, theVPars);
339 }
340
341 //=======================================================================
342 //function : FillArrayOfPnt
343 //purpose : Compute points on one surface and fill an array of points
344 //=======================================================================
FillArrayOfPnt(const Standard_Integer SurfID)345 void IntPolyh_MaillageAffinage::FillArrayOfPnt
346 (const Standard_Integer SurfID)
347 {
348 // Make sampling
349 TColStd_Array1OfReal aUpars, aVpars;
350 MakeSampling(SurfID, aUpars, aVpars);
351 // Fill array of points
352 FillArrayOfPnt(SurfID, aUpars, aVpars);
353 }
354 //=======================================================================
355 //function : FillArrayOfPnt
356 //purpose : Compute points on one surface and fill an array of points
357 // FILL AN ARRAY OF POINTS
358 //=======================================================================
FillArrayOfPnt(const Standard_Integer SurfID,const Standard_Boolean isShiftFwd)359 void IntPolyh_MaillageAffinage::FillArrayOfPnt
360 (const Standard_Integer SurfID,
361 const Standard_Boolean isShiftFwd)
362 {
363 // Make sampling
364 TColStd_Array1OfReal aUpars, aVpars;
365 MakeSampling(SurfID, aUpars, aVpars);
366 // Fill array of points
367 FillArrayOfPnt(SurfID, isShiftFwd, aUpars, aVpars);
368 }
369 //=======================================================================
370 //function : FillArrayOfPnt
371 //purpose : Compute points on one surface and fill an array of points
372 //=======================================================================
FillArrayOfPnt(const Standard_Integer SurfID,const TColStd_Array1OfReal & Upars,const TColStd_Array1OfReal & Vpars,const Standard_Real * theDeflTol)373 void IntPolyh_MaillageAffinage::FillArrayOfPnt
374 (const Standard_Integer SurfID,
375 const TColStd_Array1OfReal& Upars,
376 const TColStd_Array1OfReal& Vpars,
377 const Standard_Real *theDeflTol)
378 {
379 Standard_Boolean bDegI, bDeg;
380 Standard_Integer aNbU, aNbV, iCnt, i, j;
381 Standard_Integer aID1, aID2, aJD1, aJD2;
382 Standard_Real aTol, aU, aV, aX, aY, aZ;
383 gp_Pnt aP;
384 //
385 aNbU=(SurfID==1)? NbSamplesU1 : NbSamplesU2;
386 aNbV=(SurfID==1)? NbSamplesV1 : NbSamplesV2;
387 Bnd_Box& aBox = (SurfID==1) ? MyBox1 : MyBox2;
388 Handle(Adaptor3d_Surface)& aS=(SurfID==1)? MaSurface1:MaSurface2;
389 IntPolyh_ArrayOfPoints &TPoints=(SurfID==1)? TPoints1:TPoints2;
390 //
391 aJD1=0;
392 aJD2=0;
393 aID1=0;
394 aID2=0;
395 DegeneratedIndex(Vpars, aNbV, aS, 1, aJD1, aJD2);
396 if (!(aJD1 || aJD2)) {
397 DegeneratedIndex(Upars, aNbU, aS, 2, aID1, aID2);
398 }
399 //
400 TPoints.Init(aNbU*aNbV);
401 iCnt=0;
402 for(i=1; i<=aNbU; ++i){
403 bDegI=(aID1==i || aID2==i);
404 aU=Upars(i);
405 for(j=1; j<=aNbV; ++j){
406 aV=Vpars(j);
407 aP=aS->Value(aU, aV);
408 aP.Coord(aX, aY, aZ);
409 IntPolyh_Point& aIP=TPoints[iCnt];
410 aIP.Set(aX, aY, aZ, aU, aV);
411 //
412 bDeg=bDegI || (aJD1==j || aJD2==j);
413 if (bDeg) {
414 aIP.SetDegenerated(bDeg);
415 }
416 ++iCnt;
417 aBox.Add(aP);
418 }
419 }
420 //
421 TPoints.SetNbItems(iCnt);
422 //
423 aTol = !theDeflTol ? IntPolyh_Tools::ComputeDeflection(aS, Upars, Vpars) : *theDeflTol;
424 aTol *= 1.2;
425
426 Standard_Real a1,a2,a3,b1,b2,b3;
427 //
428 aBox.Get(a1,a2,a3,b1,b2,b3);
429 aBox.Update(a1-aTol,a2-aTol,a3-aTol,b1+aTol,b2+aTol,b3+aTol);
430 aBox.Enlarge(MyTolerance);
431 }
432
433 //=======================================================================
434 //function : FillArrayOfPnt
435 //purpose :
436 //=======================================================================
FillArrayOfPnt(const Standard_Integer SurfID,const Standard_Boolean isShiftFwd,const IntPolyh_ArrayOfPointNormal & thePointsNorm,const TColStd_Array1OfReal & theUPars,const TColStd_Array1OfReal & theVPars,const Standard_Real theDeflTol)437 void IntPolyh_MaillageAffinage::FillArrayOfPnt(const Standard_Integer SurfID,
438 const Standard_Boolean isShiftFwd,
439 const IntPolyh_ArrayOfPointNormal& thePointsNorm,
440 const TColStd_Array1OfReal& theUPars,
441 const TColStd_Array1OfReal& theVPars,
442 const Standard_Real theDeflTol)
443 {
444 Handle(Adaptor3d_Surface) aS = (SurfID == 1) ? MaSurface1 : MaSurface2;
445 IntPolyh_ArrayOfPoints& TPoints = (SurfID == 1) ? TPoints1 : TPoints2;
446 Standard_Integer aNbU = (SurfID == 1) ? NbSamplesU1 : NbSamplesU2;
447 Standard_Integer aNbV = (SurfID == 1) ? NbSamplesV1 : NbSamplesV2;
448 Bnd_Box& aBox = (SurfID==1) ? MyBox1 : MyBox2;
449
450 Standard_Integer aJD1(0), aJD2(0), aID1(0), aID2(0);
451 DegeneratedIndex(theVPars, aNbV, aS, 1, aJD1, aJD2);
452 if (!(aJD1 || aJD2))
453 DegeneratedIndex(theUPars, aNbU, aS, 2, aID1, aID2);
454
455 Standard_Boolean bDegI, bDeg;
456 Standard_Integer iCnt(0), i, j;
457 Standard_Real aX, aY, aZ, aU, aV;
458
459 TPoints.Init(thePointsNorm.NbItems());
460
461 for (i = 1; i <= aNbU; ++i)
462 {
463 aU = theUPars(i);
464 bDegI = (aID1 == i || aID2 == i);
465 for (j = 1; j <= aNbV; ++j)
466 {
467 aV = theVPars(j);
468
469 const IntPolyh_PointNormal& aPN = thePointsNorm.Value(iCnt);
470 gp_Vec aNorm = aPN.Normal.Multiplied(1.5*theDeflTol);
471 if (!isShiftFwd)
472 aNorm.Reverse();
473 gp_Pnt aP = aPN.Point.Translated(aNorm);
474
475 IntPolyh_Point& aIP = TPoints[iCnt];
476 aP.Coord(aX, aY, aZ);
477 aIP.Set(aX, aY, aZ, aU, aV);
478 bDeg = bDegI || (aJD1 == j || aJD2 == j);
479 if (bDeg)
480 aIP.SetDegenerated(bDeg);
481
482 ++iCnt;
483 aBox.Add(aP);
484 }
485 }
486
487 TPoints.SetNbItems(iCnt);
488
489 // Update box
490 Standard_Real Tol = theDeflTol*1.2;
491 Standard_Real a1,a2,a3,b1,b2,b3;
492 aBox.Get(a1,a2,a3,b1,b2,b3);
493 aBox.Update(a1-Tol,a2-Tol,a3-Tol,b1+Tol,b2+Tol,b3+Tol);
494 aBox.Enlarge(MyTolerance);
495 }
496
497 //=======================================================================
498 //function : FillArrayOfPnt
499 //purpose : Compute points on one surface and fill an array of points
500 //=======================================================================
FillArrayOfPnt(const Standard_Integer SurfID,const Standard_Boolean isShiftFwd,const TColStd_Array1OfReal & Upars,const TColStd_Array1OfReal & Vpars,const Standard_Real * theDeflTol)501 void IntPolyh_MaillageAffinage::FillArrayOfPnt
502 (const Standard_Integer SurfID,
503 const Standard_Boolean isShiftFwd,
504 const TColStd_Array1OfReal& Upars,
505 const TColStd_Array1OfReal& Vpars,
506 const Standard_Real *theDeflTol)
507 {
508 Handle(Adaptor3d_Surface) aS = (SurfID == 1) ? MaSurface1 : MaSurface2;
509 // Compute the tolerance
510 Standard_Real aTol = theDeflTol != NULL ? * theDeflTol :
511 IntPolyh_Tools::ComputeDeflection(aS, Upars, Vpars);
512 // Fill array of point normal
513 IntPolyh_ArrayOfPointNormal aPoints;
514 IntPolyh_Tools::FillArrayOfPointNormal(aS, Upars, Vpars, aPoints);
515
516 // Fill array of points
517 FillArrayOfPnt(1, isShiftFwd, aPoints, Upars, Vpars, aTol);
518 }
519
520 //=======================================================================
521 //function : CommonBox
522 //purpose :
523 //=======================================================================
CommonBox()524 void IntPolyh_MaillageAffinage::CommonBox()
525 {
526 Standard_Real XMin, YMin, ZMin, XMax, YMax, ZMax;
527 CommonBox(GetBox(1), GetBox(2), XMin, YMin, ZMin, XMax, YMax, ZMax);
528 }
529
530 //=======================================================================
531 //function : CommonBox
532 //purpose : Compute the common box witch is the intersection
533 // of the two bounding boxes, and mark the points of
534 // the two surfaces that are inside.
535 // REJECTION BOUNDING BOXES
536 // DETERMINATION OF THE COMMON BOX
537 //=======================================================================
CommonBox(const Bnd_Box &,const Bnd_Box &,Standard_Real & XMin,Standard_Real & YMin,Standard_Real & ZMin,Standard_Real & XMax,Standard_Real & YMax,Standard_Real & ZMax)538 void IntPolyh_MaillageAffinage::CommonBox (const Bnd_Box &,
539 const Bnd_Box &,
540 Standard_Real &XMin,
541 Standard_Real &YMin,
542 Standard_Real &ZMin,
543 Standard_Real &XMax,
544 Standard_Real &YMax,
545 Standard_Real &ZMax)
546 {
547 Standard_Real x10,y10,z10,x11,y11,z11;
548 Standard_Real x20,y20,z20,x21,y21,z21;
549
550 MyBox1.Get(x10,y10,z10,x11,y11,z11);
551 MyBox2.Get(x20,y20,z20,x21,y21,z21);
552 XMin = 0.;
553 YMin = 0.;
554 ZMin = 0.;
555 XMax = 0.;
556 YMax = 0.;
557 ZMax = 0.;
558
559 if((x10>x21)||(x20>x11)||(y10>y21)||(y20>y11)||(z10>z21)||(z20>z11)) {
560 }
561 else {
562 if(x11<=x21) XMax=x11; else { if(x21<=x11) XMax=x21;}
563 if(x20<=x10) XMin=x10; else { if(x10<=x20) XMin=x20;}
564 if(y11<=y21) YMax=y11; else { if(y21<=y11) YMax=y21;}
565 if(y20<=y10) YMin=y10; else { if(y10<=y20) YMin=y20;}
566 if(z11<=z21) ZMax=z11; else { if(z21<=z11) ZMax=z21;}
567 if(z20<=z10) ZMin=z10; else { if(z10<=z20) ZMin=z20;}
568
569 if(((XMin==XMax)&&(!(YMin==YMax)&&!(ZMin==ZMax)))
570 ||((YMin==YMax)&&(!(XMin==XMax)&&!(ZMin==ZMax)))//ou exclusif ??
571 ||((ZMin==ZMax)&&(!(XMin==XMax)&&!(YMin==YMax)))) {
572 }
573 }
574 //
575 Standard_Real X,Y,Z;
576 X=XMax-XMin;
577 Y=YMax-YMin;
578 Z=ZMax-ZMin;
579 //extension of the box
580 if( (X==0)&&(Y!=0) ) X=Y*0.1;
581 else if( (X==0)&&(Z!=0) ) X=Z*0.1;
582 else X*=0.1;
583
584 if( (Y==0)&&(X!=0) ) Y=X*0.1;
585 else if( (Y==0)&&(Z!=0) ) Y=Z*0.1;
586 else Y*=0.1;
587
588 if( (Z==0)&&(X!=0) ) Z=X*0.1;
589 else if( (Z==0)&&(Y!=0) ) Z=Y*0.1;
590 else Z*=0.1;
591
592
593 if( (X==0)&&(Y==0)&&(Z==0) ) {
594
595
596 }
597 XMin-=X; XMax+=X;
598 YMin-=Y; YMax+=Y;
599 ZMin-=Z; ZMax+=Z;
600
601 //Marking of points included in the common
602 const Standard_Integer FinTP1 = TPoints1.NbItems();
603 // for(Standard_Integer i=0; i<FinTP1; i++) {
604 Standard_Integer i ;
605 for( i=0; i<FinTP1; i++) {
606 IntPolyh_Point & Pt1 = TPoints1[i];
607 Standard_Integer r;
608 if(Pt1.X()<XMin) {
609 r=1;
610 }
611 else {
612 if(Pt1.X()>XMax) {
613 r=2;
614 } else {
615 r=0;
616 }
617 }
618 if(Pt1.Y()<YMin) {
619 r|=4;
620 }
621 else {
622 if(Pt1.Y()>YMax) {
623 r|=8;
624 }
625 }
626 if(Pt1.Z()<ZMin) {
627 r|=16;
628 } else {
629 if(Pt1.Z()>ZMax) {
630 r|=32;
631 }
632 }
633 Pt1.SetPartOfCommon(r);
634 }
635
636 const Standard_Integer FinTP2 = TPoints2.NbItems();
637 for(Standard_Integer ii=0; ii<FinTP2; ii++) {
638 IntPolyh_Point & Pt2 = TPoints2[ii];
639 Standard_Integer rr;
640 if(Pt2.X()<XMin) {
641 rr=1;
642 }
643 else {
644 if(Pt2.X()>XMax) {
645 rr=2;
646 } else {
647 rr=0;
648 }
649 }
650 if(Pt2.Y()<YMin) {
651 rr|=4;
652 }
653 else {
654 if(Pt2.Y()>YMax) {
655 rr|=8;
656 }
657 }
658 if(Pt2.Z()<ZMin) {
659 rr|=16;
660 }
661 else {
662 if(Pt2.Z()>ZMax) {
663 rr|=32;
664 }
665 }
666 Pt2.SetPartOfCommon(rr);
667 }
668 }
669 //=======================================================================
670 //function : FillArrayOfEdges
671 //purpose : Compute edges from the array of points
672 // FILL THE ARRAY OF EDGES
673 //=======================================================================
FillArrayOfEdges(const Standard_Integer SurfID)674 void IntPolyh_MaillageAffinage::FillArrayOfEdges
675 (const Standard_Integer SurfID)
676 {
677
678 IntPolyh_ArrayOfEdges &TEdges=(SurfID==1)? TEdges1:TEdges2;
679 IntPolyh_ArrayOfTriangles &TTriangles=(SurfID==1)? TTriangles1:TTriangles2;
680 Standard_Integer NbSamplesU=(SurfID==1)? NbSamplesU1:NbSamplesU2;
681 Standard_Integer NbSamplesV=(SurfID==1)? NbSamplesV1:NbSamplesV2;
682
683 //NbEdges = 3 + 3*(NbSamplesV-2) + 3*(NbSamplesU-2) +
684 // + 3*(NbSamplesU-2)*(NbSamplesV-2) + (NbSamplesV-1) + (NbSamplesU-1);
685 //NbSamplesU and NbSamples cannot be less than 2, so
686 Standard_Integer NbEdges = 3*NbSamplesU*NbSamplesV - 2*(NbSamplesU+NbSamplesV) + 1;
687 TEdges.Init(NbEdges);
688
689 Standard_Integer CpteurTabEdges=0;
690
691 //maillage u0 v0
692 TEdges[CpteurTabEdges].SetFirstPoint(0); // U V
693 TEdges[CpteurTabEdges].SetSecondPoint(1); // U V+1
694 TEdges[CpteurTabEdges].SetSecondTriangle(0);
695 TTriangles[0].SetEdgeAndOrientation(TEdges[CpteurTabEdges], CpteurTabEdges);
696 CpteurTabEdges++;
697
698 TEdges[CpteurTabEdges].SetFirstPoint(0); // U V
699 TEdges[CpteurTabEdges].SetSecondPoint(NbSamplesV); // U+1 V
700 TEdges[CpteurTabEdges].SetFirstTriangle(1);
701 TTriangles[1].SetEdgeAndOrientation(TEdges[CpteurTabEdges], CpteurTabEdges);
702 CpteurTabEdges++;
703
704 TEdges[CpteurTabEdges].SetFirstPoint(0); // U V
705 TEdges[CpteurTabEdges].SetSecondPoint(NbSamplesV+1); // U+1 V+1
706 TEdges[CpteurTabEdges].SetFirstTriangle(0);
707 TTriangles[0].SetEdgeAndOrientation(TEdges[CpteurTabEdges], CpteurTabEdges);
708 TEdges[CpteurTabEdges].SetSecondTriangle(1);
709 TTriangles[1].SetEdgeAndOrientation(TEdges[CpteurTabEdges], CpteurTabEdges);
710 CpteurTabEdges++;
711
712 //maillage surU=u0
713 Standard_Integer PntInit=1;
714 Standard_Integer BoucleMeshV;
715 for(BoucleMeshV=1; BoucleMeshV<NbSamplesV-1;BoucleMeshV++){
716 TEdges[CpteurTabEdges].SetFirstPoint(PntInit); // U V
717 TEdges[CpteurTabEdges].SetSecondPoint(PntInit+1); // U V+1
718 // TEdges[CpteurTabEdges].SetFirstTriangle(-1);
719 TEdges[CpteurTabEdges].SetSecondTriangle(BoucleMeshV*2);
720 TTriangles[BoucleMeshV*2].SetEdgeAndOrientation(TEdges[CpteurTabEdges], CpteurTabEdges);
721 CpteurTabEdges++;
722
723 TEdges[CpteurTabEdges].SetFirstPoint(PntInit); // U V
724 TEdges[CpteurTabEdges].SetSecondPoint(PntInit+NbSamplesV+1); // U+1 V+1
725 TEdges[CpteurTabEdges].SetFirstTriangle(BoucleMeshV*2);
726 TTriangles[BoucleMeshV*2].SetEdgeAndOrientation(TEdges[CpteurTabEdges], CpteurTabEdges);
727 TEdges[CpteurTabEdges].SetSecondTriangle(BoucleMeshV*2+1);
728 TTriangles[BoucleMeshV*2+1].SetEdgeAndOrientation(TEdges[CpteurTabEdges], CpteurTabEdges);
729 CpteurTabEdges++;
730
731 TEdges[CpteurTabEdges].SetFirstPoint(PntInit); // U V
732 TEdges[CpteurTabEdges].SetSecondPoint(PntInit+NbSamplesV); // U+1 V
733 TEdges[CpteurTabEdges].SetFirstTriangle(BoucleMeshV*2+1);
734 TTriangles[BoucleMeshV*2+1].SetEdgeAndOrientation(TEdges[CpteurTabEdges], CpteurTabEdges);
735 TEdges[CpteurTabEdges].SetSecondTriangle(BoucleMeshV*2-2);
736 TTriangles[BoucleMeshV*2-2].SetEdgeAndOrientation(TEdges[CpteurTabEdges], CpteurTabEdges);
737 CpteurTabEdges++;
738 PntInit++;
739 }
740
741 //maillage sur V=v0
742 PntInit=NbSamplesV;
743 for(BoucleMeshV=1; BoucleMeshV<NbSamplesU-1;BoucleMeshV++){
744 TEdges[CpteurTabEdges].SetFirstPoint(PntInit); // U V
745 TEdges[CpteurTabEdges].SetSecondPoint(PntInit+1); // U V+1
746 TEdges[CpteurTabEdges].SetFirstTriangle((BoucleMeshV-1)*(NbSamplesV-1)*2+1);
747 TTriangles[(BoucleMeshV-1)*(NbSamplesV-1)*2+1].SetEdgeAndOrientation(TEdges[CpteurTabEdges], CpteurTabEdges);
748 TEdges[CpteurTabEdges].SetSecondTriangle(BoucleMeshV*(NbSamplesV-1)*2);
749 TTriangles[BoucleMeshV*(NbSamplesV-1)*2].SetEdgeAndOrientation(TEdges[CpteurTabEdges], CpteurTabEdges);
750 CpteurTabEdges++;
751
752 TEdges[CpteurTabEdges].SetFirstPoint(PntInit); // U V
753 TEdges[CpteurTabEdges].SetSecondPoint(PntInit+NbSamplesV+1); // U+1 V+1
754 TEdges[CpteurTabEdges].SetFirstTriangle(BoucleMeshV*(NbSamplesV-1)*2);
755 TTriangles[BoucleMeshV*(NbSamplesV-1)*2].SetEdgeAndOrientation(TEdges[CpteurTabEdges], CpteurTabEdges);
756 TEdges[CpteurTabEdges].SetSecondTriangle(BoucleMeshV*(NbSamplesV-1)*2+1);
757 TTriangles[BoucleMeshV*(NbSamplesV-1)*2+1].SetEdgeAndOrientation(TEdges[CpteurTabEdges], CpteurTabEdges);
758 CpteurTabEdges++;
759
760 TEdges[CpteurTabEdges].SetFirstPoint(PntInit); // U V
761 TEdges[CpteurTabEdges].SetSecondPoint(PntInit+NbSamplesV); // U+1 V
762 TEdges[CpteurTabEdges].SetFirstTriangle(BoucleMeshV*(NbSamplesV-1)*2+1);
763 TTriangles[BoucleMeshV*(NbSamplesV-1)*2+1].SetEdgeAndOrientation(TEdges[CpteurTabEdges], CpteurTabEdges);
764 CpteurTabEdges++;
765 PntInit+=NbSamplesV;
766 }
767
768 PntInit=NbSamplesV+1;
769 //To provide recursion I associate a point with three edges
770 for(Standard_Integer BoucleMeshU=1; BoucleMeshU<NbSamplesU-1; BoucleMeshU++){
771 for(BoucleMeshV=1; BoucleMeshV<NbSamplesV-1;BoucleMeshV++){
772 TEdges[CpteurTabEdges].SetFirstPoint(PntInit); // U V
773 TEdges[CpteurTabEdges].SetSecondPoint(PntInit+1); // U V+1
774 TEdges[CpteurTabEdges].SetFirstTriangle((NbSamplesV-1)*2*(BoucleMeshU-1)+BoucleMeshV*2+1);
775 TTriangles[(NbSamplesV-1)*2*(BoucleMeshU-1)+BoucleMeshV*2+1].SetEdgeAndOrientation(TEdges[CpteurTabEdges], CpteurTabEdges);
776 TEdges[CpteurTabEdges].SetSecondTriangle((NbSamplesV-1)*2*BoucleMeshU+BoucleMeshV*2);
777 TTriangles[(NbSamplesV-1)*2*BoucleMeshU+BoucleMeshV*2].SetEdgeAndOrientation(TEdges[CpteurTabEdges], CpteurTabEdges);
778 CpteurTabEdges++;
779
780 TEdges[CpteurTabEdges].SetFirstPoint(PntInit); // U V
781 TEdges[CpteurTabEdges].SetSecondPoint(PntInit+NbSamplesV+1); // U+1 V+1
782 TEdges[CpteurTabEdges].SetFirstTriangle((NbSamplesV-1)*2*BoucleMeshU+BoucleMeshV*2);
783 TTriangles[(NbSamplesV-1)*2*BoucleMeshU+BoucleMeshV*2].SetEdgeAndOrientation(TEdges[CpteurTabEdges], CpteurTabEdges);
784 TEdges[CpteurTabEdges].SetSecondTriangle((NbSamplesV-1)*2*BoucleMeshU+BoucleMeshV*2+1);
785 TTriangles[(NbSamplesV-1)*2*BoucleMeshU+BoucleMeshV*2+1].SetEdgeAndOrientation(TEdges[CpteurTabEdges], CpteurTabEdges);
786 CpteurTabEdges++;
787
788 TEdges[CpteurTabEdges].SetFirstPoint(PntInit); // U V
789 TEdges[CpteurTabEdges].SetSecondPoint(PntInit+NbSamplesV); // U+1 V
790 TEdges[CpteurTabEdges].SetFirstTriangle((NbSamplesV-1)*2*BoucleMeshU+BoucleMeshV*2+1);
791 TTriangles[(NbSamplesV-1)*2*BoucleMeshU+BoucleMeshV*2+1].SetEdgeAndOrientation(TEdges[CpteurTabEdges], CpteurTabEdges);
792 TEdges[CpteurTabEdges].SetSecondTriangle((NbSamplesV-1)*2*BoucleMeshU+BoucleMeshV*2-2);
793 TTriangles[(NbSamplesV-1)*2*BoucleMeshU+BoucleMeshV*2-2].SetEdgeAndOrientation(TEdges[CpteurTabEdges], CpteurTabEdges);
794 CpteurTabEdges++;
795 PntInit++;//Pass to the next point
796 }
797 PntInit++;//Pass the last point of the column
798 PntInit++;//Pass the first point of the next column
799 }
800
801 //close mesh on U=u1
802 PntInit=(NbSamplesU-1)*NbSamplesV; //point U=u1 V=0
803 for(BoucleMeshV=0; BoucleMeshV<NbSamplesV-1; BoucleMeshV++){
804 TEdges[CpteurTabEdges].SetFirstPoint(PntInit); //U=u1 V
805 TEdges[CpteurTabEdges].SetSecondPoint(PntInit+1); //U=u1 V+1
806 TEdges[CpteurTabEdges].SetFirstTriangle((NbSamplesU-2)*(NbSamplesV-1)*2+BoucleMeshV*2+1);
807 TTriangles[(NbSamplesU-2)*(NbSamplesV-1)*2+BoucleMeshV*2+1].SetEdgeAndOrientation(TEdges[CpteurTabEdges], CpteurTabEdges);
808 CpteurTabEdges++;
809 PntInit++;
810 }
811
812 //close mesh on V=v1
813 for(BoucleMeshV=0; BoucleMeshV<NbSamplesU-1;BoucleMeshV++){
814 TEdges[CpteurTabEdges].SetFirstPoint(NbSamplesV-1+BoucleMeshV*NbSamplesV); // U V=v1
815 TEdges[CpteurTabEdges].SetSecondPoint(NbSamplesV-1+(BoucleMeshV+1)*NbSamplesV); //U+1 V=v1
816 TEdges[CpteurTabEdges].SetSecondTriangle(BoucleMeshV*2*(NbSamplesV-1)+(NbSamplesV-2)*2);
817 TTriangles[BoucleMeshV*2*(NbSamplesV-1)+(NbSamplesV-2)*2].SetEdgeAndOrientation(TEdges[CpteurTabEdges], CpteurTabEdges);
818 CpteurTabEdges++;
819 }
820 TEdges.SetNbItems(CpteurTabEdges);
821 }
822
823 //=======================================================================
824 //function : FillArrayOfTriangles
825 //purpose : Compute triangles from the array of points, and --
826 // mark the triangles that use marked points by the
827 // CommonBox function.
828 // FILL THE ARRAY OF TRIANGLES
829 //=======================================================================
FillArrayOfTriangles(const Standard_Integer SurfID)830 void IntPolyh_MaillageAffinage::FillArrayOfTriangles
831 (const Standard_Integer SurfID)
832 {
833 Standard_Integer CpteurTabTriangles=0;
834 Standard_Integer PntInit=0;
835
836 IntPolyh_ArrayOfPoints &TPoints=(SurfID==1)? TPoints1:TPoints2;
837 IntPolyh_ArrayOfTriangles &TTriangles=(SurfID==1)? TTriangles1:TTriangles2;
838 Standard_Integer NbSamplesU=(SurfID==1)? NbSamplesU1:NbSamplesU2;
839 Standard_Integer NbSamplesV=(SurfID==1)? NbSamplesV1:NbSamplesV2;
840
841 TTriangles.Init(2*(NbSamplesU-1)*(NbSamplesV-1));
842 //To provide recursion, I associate a point with two triangles
843 for(Standard_Integer BoucleMeshU=0; BoucleMeshU<NbSamplesU-1; BoucleMeshU++){
844 for(Standard_Integer BoucleMeshV=0; BoucleMeshV<NbSamplesV-1;BoucleMeshV++){
845
846 // FIRST TRIANGLE
847 TTriangles[CpteurTabTriangles].SetFirstPoint(PntInit); // U V
848 TTriangles[CpteurTabTriangles].SetSecondPoint(PntInit+1); // U V+1
849 TTriangles[CpteurTabTriangles].SetThirdPoint(PntInit+NbSamplesV+1); // U+1 V+1
850
851 // IF ITS EDGE CONTACTS WITH THE COMMON BOX IP REMAINS = A 1
852 if( ( (TPoints[PntInit].PartOfCommon()) & (TPoints[PntInit+1].PartOfCommon()) )
853 &&( (TPoints[PntInit+1].PartOfCommon()) & (TPoints[PntInit+NbSamplesV+1].PartOfCommon()))
854 &&( (TPoints[PntInit+NbSamplesV+1].PartOfCommon()) & (TPoints[PntInit].PartOfCommon())) )
855 //IF NOT IP=0
856 TTriangles[CpteurTabTriangles].SetIntersectionPossible(Standard_False);
857
858 CpteurTabTriangles++;
859
860 //SECOND TRIANGLE
861 TTriangles[CpteurTabTriangles].SetFirstPoint(PntInit); // U V
862 TTriangles[CpteurTabTriangles].SetSecondPoint(PntInit+NbSamplesV+1); // U+1 V+1
863 TTriangles[CpteurTabTriangles].SetThirdPoint(PntInit+NbSamplesV); // U+1 V
864
865
866 if( ( (TPoints[PntInit].PartOfCommon()) & (TPoints[PntInit+NbSamplesV+1].PartOfCommon()) )
867 &&( (TPoints[PntInit+NbSamplesV+1].PartOfCommon()) & (TPoints[PntInit+NbSamplesV].PartOfCommon()))
868 &&( (TPoints[PntInit+NbSamplesV].PartOfCommon()) & (TPoints[PntInit].PartOfCommon())) )
869 TTriangles[CpteurTabTriangles].SetIntersectionPossible(Standard_False);
870
871
872 CpteurTabTriangles++;
873
874 PntInit++;//Pass to the next point
875 }
876 PntInit++;//Pass the last point of the column
877 }
878 TTriangles.SetNbItems(CpteurTabTriangles);
879 const Standard_Integer FinTT = TTriangles.NbItems();
880 if (FinTT==0) {
881 }
882 }
883 //=======================================================================
884 //function : CommonPartRefinement
885 //purpose : Refine systematicaly all marked triangles of both surfaces
886 // REFINING OF THE COMMON
887 //=======================================================================
CommonPartRefinement()888 void IntPolyh_MaillageAffinage::CommonPartRefinement()
889 {
890 Standard_Integer FinInit1 = TTriangles1.NbItems();
891 for(Standard_Integer i=0; i<FinInit1; i++) {
892 if(TTriangles1[i].IsIntersectionPossible())
893 TTriangles1[i].MiddleRefinement(i,MaSurface1,TPoints1,TTriangles1,TEdges1);
894 }
895
896 Standard_Integer FinInit2=TTriangles2.NbItems();
897 for(Standard_Integer ii=0; ii<FinInit2; ii++) {
898 if(TTriangles2[ii].IsIntersectionPossible())
899 TTriangles2[ii].MiddleRefinement(ii,MaSurface2,TPoints2,TTriangles2,TEdges2);
900 }
901
902 }
903 //=======================================================================
904 //function : LocalSurfaceRefinement
905 //purpose : Refine systematicaly all marked triangles of ONE surface
906 //=======================================================================
LocalSurfaceRefinement(const Standard_Integer SurfID)907 void IntPolyh_MaillageAffinage::LocalSurfaceRefinement(const Standard_Integer SurfID) {
908 //refine locally, but systematically the chosen surface
909 if (SurfID==1) {
910 const Standard_Integer FinInit1 = TTriangles1.NbItems();
911 for(Standard_Integer i=0; i<FinInit1; i++) {
912 if(TTriangles1[i].IsIntersectionPossible())
913 TTriangles1[i].MiddleRefinement(i,MaSurface1,TPoints1,TTriangles1,TEdges1);
914 }
915 }
916 //
917 if (SurfID==2) {
918 const Standard_Integer FinInit2 = TTriangles2.NbItems();
919 for(Standard_Integer ii=0; ii<FinInit2; ii++) {
920 if(TTriangles2[ii].IsIntersectionPossible())
921 TTriangles2[ii].MiddleRefinement(ii,MaSurface2,TPoints2,TTriangles2,TEdges2);
922 }
923 }
924 }
925 //=======================================================================
926 //function : ComputeDeflections
927 //purpose : Compute deflection for all triangles of one
928 // surface,and sort min and max of deflections
929 // REFINING PART
930 // Calculation of the deflection of all triangles
931 // --> deflection max
932 // --> deflection min
933 //=======================================================================
ComputeDeflections(const Standard_Integer SurfID)934 void IntPolyh_MaillageAffinage::ComputeDeflections
935 (const Standard_Integer SurfID)
936 {
937 Handle(Adaptor3d_Surface) aSurface=(SurfID==1)? MaSurface1:MaSurface2;
938 IntPolyh_ArrayOfPoints &TPoints=(SurfID==1)? TPoints1:TPoints2;
939 IntPolyh_ArrayOfTriangles &TTriangles=(SurfID==1)? TTriangles1:TTriangles2;
940 Standard_Real &FlecheMin=(SurfID==1)? FlecheMin1:FlecheMin2;
941 Standard_Real &FlecheMax=(SurfID==1)? FlecheMax1:FlecheMax2;
942
943 FlecheMax=-RealLast();
944 FlecheMin=RealLast();
945 const Standard_Integer FinTT = TTriangles.NbItems();
946
947 for(Standard_Integer i = 0; i < FinTT; i++) {
948 IntPolyh_Triangle& aTriangle = TTriangles[i];
949 Standard_Real Fleche = aTriangle.ComputeDeflection(aSurface, TPoints);
950 if (Fleche > FlecheMax)
951 FlecheMax = Fleche;
952 if (Fleche < FlecheMin)
953 FlecheMin = Fleche;
954 }
955 }
956
957 //=======================================================================
958 //function : TrianglesDeflectionsRefinement
959 //purpose : Refinement of the triangles depending on the deflection
960 //=======================================================================
961 static
TrianglesDeflectionsRefinement(const Handle (Adaptor3d_Surface)& theS1,IntPolyh_ArrayOfTriangles & theTriangles1,IntPolyh_ArrayOfEdges & theEdges1,IntPolyh_ArrayOfPoints & thePoints1,const Standard_Real theFlecheCritique1,const Handle (Adaptor3d_Surface)& theS2,IntPolyh_ArrayOfTriangles & theTriangles2,IntPolyh_ArrayOfEdges & theEdges2,IntPolyh_ArrayOfPoints & thePoints2,const Standard_Real theFlecheCritique2)962 void TrianglesDeflectionsRefinement(const Handle(Adaptor3d_Surface)& theS1,
963 IntPolyh_ArrayOfTriangles& theTriangles1,
964 IntPolyh_ArrayOfEdges& theEdges1,
965 IntPolyh_ArrayOfPoints& thePoints1,
966 const Standard_Real theFlecheCritique1,
967 const Handle(Adaptor3d_Surface)& theS2,
968 IntPolyh_ArrayOfTriangles& theTriangles2,
969 IntPolyh_ArrayOfEdges& theEdges2,
970 IntPolyh_ArrayOfPoints& thePoints2,
971 const Standard_Real theFlecheCritique2)
972 {
973 // Find the intersecting triangles
974 IntPolyh_IndexedDataMapOfIntegerListOfInteger aDMILI;
975 GetInterferingTriangles(theTriangles1, thePoints1, theTriangles2, thePoints2, aDMILI);
976 //
977 // Interfering triangles of second surface
978 TColStd_MapOfInteger aMIntS2;
979 // Since the number of the triangles may change during refinement,
980 // save the number of triangles before refinement
981 Standard_Integer FinTT1 = theTriangles1.NbItems();
982 Standard_Integer FinTT2 = theTriangles2.NbItems();
983 //
984 // Analyze interfering triangles
985 for (Standard_Integer i_S1 = 0; i_S1 < FinTT1; i_S1++) {
986 IntPolyh_Triangle& aTriangle1 = theTriangles1[i_S1];
987 if (!aTriangle1.IsIntersectionPossible()) {
988 continue;
989 }
990 //
991 const TColStd_ListOfInteger *pLI = aDMILI.Seek(i_S1);
992 if (!pLI || pLI->IsEmpty()) {
993 // Mark non-interfering triangles of S1 to avoid their repeated usage
994 aTriangle1.SetIntersectionPossible(Standard_False);
995 continue;
996 }
997 //
998 if (aTriangle1.Deflection() > theFlecheCritique1) {
999 aTriangle1.MiddleRefinement(i_S1, theS1, thePoints1, theTriangles1, theEdges1);
1000 }
1001 //
1002 TColStd_ListOfInteger::Iterator Iter(*pLI);
1003 for (; Iter.More(); Iter.Next()) {
1004 Standard_Integer i_S2 = Iter.Value();
1005 if (aMIntS2.Add(i_S2)) {
1006 IntPolyh_Triangle & aTriangle2 = theTriangles2[i_S2];
1007 if (aTriangle2.Deflection() > theFlecheCritique2) {
1008 // Refinement of the larger triangles
1009 aTriangle2.MiddleRefinement(i_S2, theS2, thePoints2, theTriangles2, theEdges2);
1010 }
1011 }
1012 }
1013 }
1014 //
1015 // Mark non-interfering triangles of S2 to avoid their repeated usage
1016 if (aMIntS2.Extent() < FinTT2) {
1017 for (Standard_Integer i_S2 = 0; i_S2 < FinTT2; i_S2++) {
1018 if (!aMIntS2.Contains(i_S2)) {
1019 theTriangles2[i_S2].SetIntersectionPossible(Standard_False);
1020 }
1021 }
1022 }
1023 }
1024 //=======================================================================
1025 //function : LargeTrianglesDeflectionsRefinement
1026 //purpose : Refinement of the large triangles in case one surface is
1027 // much smaller then the other.
1028 //=======================================================================
1029 static
LargeTrianglesDeflectionsRefinement(const Handle (Adaptor3d_Surface)& theS,IntPolyh_ArrayOfTriangles & theTriangles,IntPolyh_ArrayOfEdges & theEdges,IntPolyh_ArrayOfPoints & thePoints,const Bnd_Box & theOppositeBox)1030 void LargeTrianglesDeflectionsRefinement(const Handle(Adaptor3d_Surface)& theS,
1031 IntPolyh_ArrayOfTriangles& theTriangles,
1032 IntPolyh_ArrayOfEdges& theEdges,
1033 IntPolyh_ArrayOfPoints& thePoints,
1034 const Bnd_Box& theOppositeBox)
1035 {
1036 // Find all triangles of the bigger surface with bounding boxes
1037 // overlapping the bounding box the other surface
1038 TColStd_ListOfInteger aLIT;
1039 Standard_Integer i, aNbT = theTriangles.NbItems();
1040 for (i = 0; i < aNbT; ++i) {
1041 IntPolyh_Triangle& aTriangle = theTriangles[i];
1042 if (!aTriangle.IsIntersectionPossible() || aTriangle.IsDegenerated()) {
1043 continue;
1044 }
1045 //
1046 const Bnd_Box& aBox = aTriangle.BoundingBox(thePoints);
1047 if (aBox.IsVoid()) {
1048 continue;
1049 }
1050 //
1051 if (aBox.IsOut(theOppositeBox)) {
1052 aTriangle.SetIntersectionPossible(Standard_False);
1053 continue;
1054 }
1055 //
1056 aLIT.Append(i);
1057 }
1058 //
1059 if (aLIT.IsEmpty()) {
1060 return;
1061 }
1062 //
1063 // One surface is very small relatively to the other.
1064 // The criterion of refining for large surface depends on the size of
1065 // the bounding box of the other - since the criterion should be minimized,
1066 // the smallest side of the bounding box is taken
1067 Standard_Real x0, y0, z0, x1, y1, z1;
1068 theOppositeBox.Get(x0, y0, z0, x1, y1, z1);
1069 Standard_Real dx = Abs(x1 - x0);
1070 Standard_Real dy = Abs(y1 - y0);
1071 Standard_Real diag = Abs(z1 - z0);
1072 Standard_Real dd = (dx > dy) ? dy : dx;
1073 if (diag > dd) diag = dd;
1074
1075 // calculation of the criterion of refining
1076 Standard_Real CritereAffinage = 0.0;
1077 Standard_Real DiagPonderation = 0.5;
1078 CritereAffinage = diag*DiagPonderation;
1079 //
1080 // The deflection of the greater is compared to the size of the smaller
1081 TColStd_ListIteratorOfListOfInteger Iter(aLIT);
1082 for (; Iter.More(); Iter.Next()) {
1083 i = Iter.Value();
1084 IntPolyh_Triangle & aTriangle = theTriangles[i];
1085 if (aTriangle.Deflection() > CritereAffinage) {
1086 aTriangle.MultipleMiddleRefinement(CritereAffinage, theOppositeBox, i,
1087 theS, thePoints, theTriangles, theEdges);
1088 }
1089 else {
1090 aTriangle.MiddleRefinement(i, theS, thePoints, theTriangles, theEdges);
1091 }
1092 }
1093 }
1094 //=======================================================================
1095 //function : TrianglesDeflectionsRefinementBSB
1096 //purpose : Refine both surfaces using UBTree as rejection.
1097 // The criterion used to refine a triangle are:
1098 // - The deflection;
1099 // - The size of the - bounding boxes (one surface may be
1100 // very small compared to the other).
1101 //=======================================================================
TrianglesDeflectionsRefinementBSB()1102 void IntPolyh_MaillageAffinage::TrianglesDeflectionsRefinementBSB()
1103 {
1104 // To estimate a surface in general it can be interesting
1105 // to calculate all deflections
1106 ComputeDeflections(1);
1107 // Check deflection at output
1108 Standard_Real FlecheCritique1;
1109 if (FlecheMin1 > FlecheMax1) {
1110 return;
1111 }
1112 else {
1113 // fleche min + (flechemax-flechemin) * 80/100
1114 FlecheCritique1 = FlecheMin1*0.2 + FlecheMax1*0.8;
1115 }
1116
1117 // Compute deflections for the second surface
1118 ComputeDeflections(2);
1119
1120 //-- Check arrows at output
1121 Standard_Real FlecheCritique2;
1122 if (FlecheMin2 > FlecheMax2) {
1123 return;
1124 }
1125 else {
1126 //fleche min + (flechemax-flechemin) * 80/100
1127 FlecheCritique2 = FlecheMin2*0.2 + FlecheMax2*0.8;
1128 }
1129
1130 // The greatest of two bounding boxes created in FillArrayOfPoints is found.
1131 // Then this value is weighted depending on the discretization
1132 // (NbSamplesU and NbSamplesV)
1133 Standard_Real diag1, diag2;
1134 Standard_Real x0, y0, z0, x1, y1, z1;
1135
1136 MyBox1.Get(x0, y0, z0, x1, y1, z1);
1137 x0 -= x1; y0 -= y1; z0 -= z1;
1138 diag1 = x0*x0 + y0*y0 + z0*z0;
1139 const Standard_Real NbSamplesUV1 = Standard_Real(NbSamplesU1) * Standard_Real(NbSamplesV1);
1140 diag1 /= NbSamplesUV1;
1141
1142 MyBox2.Get(x0, y0, z0, x1, y1, z1);
1143 x0 -= x1; y0 -= y1; z0 -= z1;
1144 diag2 = x0*x0 + y0*y0 + z0*z0;
1145 const Standard_Real NbSamplesUV2 = Standard_Real(NbSamplesU2) * Standard_Real(NbSamplesV2);
1146 diag2 /= NbSamplesUV2;
1147
1148 // The surface with the greatest bounding box is "discretized"
1149 if (diag1 < diag2) {
1150 // second is discretized
1151 if (FlecheCritique2 < diag1) {
1152 // The corresponding sizes are not too disproportional
1153 TrianglesDeflectionsRefinement(MaSurface1, TTriangles1, TEdges1, TPoints1, FlecheCritique1,
1154 MaSurface2, TTriangles2, TEdges2, TPoints2, FlecheCritique2);
1155 }
1156 else {
1157 // second surface is much larger then the first
1158 LargeTrianglesDeflectionsRefinement(MaSurface2, TTriangles2, TEdges2, TPoints2, MyBox1);
1159 }
1160 }
1161 else {
1162 // first is discretized
1163 if (FlecheCritique1 < diag2) {
1164 // The corresponding sizes are not too disproportional
1165 TrianglesDeflectionsRefinement(MaSurface2, TTriangles2, TEdges2, TPoints2, FlecheCritique2,
1166 MaSurface1, TTriangles1, TEdges1, TPoints1, FlecheCritique1);
1167 }
1168 else {
1169 // first surface is much larger then the second
1170 LargeTrianglesDeflectionsRefinement(MaSurface1, TTriangles1, TEdges1, TPoints1, MyBox2);
1171 }
1172 }
1173 }
1174 //=======================================================================
1175 //function : maxSR
1176 //purpose : This function is used for the function project6
1177 //=======================================================================
maxSR(const Standard_Real a,const Standard_Real b,const Standard_Real c)1178 inline Standard_Real maxSR(const Standard_Real a,
1179 const Standard_Real b,
1180 const Standard_Real c)
1181 {
1182 Standard_Real t = a;
1183 if (b > t) t = b;
1184 if (c > t) t = c;
1185 return t;
1186 }
1187 //=======================================================================
1188 //function : minSR
1189 //purpose : This function is used for the function project6
1190 //=======================================================================
minSR(const Standard_Real a,const Standard_Real b,const Standard_Real c)1191 inline Standard_Real minSR(const Standard_Real a,
1192 const Standard_Real b,
1193 const Standard_Real c)
1194 {
1195 Standard_Real t = a;
1196 if (b < t) t = b;
1197 if (c < t) t = c;
1198 return t;
1199 }
1200 //=======================================================================
1201 //function : project6
1202 //purpose : This function is used for the function TriContact
1203 //=======================================================================
project6(const IntPolyh_Point & ax,const IntPolyh_Point & p1,const IntPolyh_Point & p2,const IntPolyh_Point & p3,const IntPolyh_Point & q1,const IntPolyh_Point & q2,const IntPolyh_Point & q3)1204 Standard_Integer project6(const IntPolyh_Point &ax,
1205 const IntPolyh_Point &p1,
1206 const IntPolyh_Point &p2,
1207 const IntPolyh_Point &p3,
1208 const IntPolyh_Point &q1,
1209 const IntPolyh_Point &q2,
1210 const IntPolyh_Point &q3)
1211 {
1212 Standard_Real P1 = ax.Dot(p1);
1213 Standard_Real P2 = ax.Dot(p2);
1214 Standard_Real P3 = ax.Dot(p3);
1215 Standard_Real Q1 = ax.Dot(q1);
1216 Standard_Real Q2 = ax.Dot(q2);
1217 Standard_Real Q3 = ax.Dot(q3);
1218
1219 Standard_Real mx1 = maxSR(P1, P2, P3);
1220 Standard_Real mn1 = minSR(P1, P2, P3);
1221 Standard_Real mx2 = maxSR(Q1, Q2, Q3);
1222 Standard_Real mn2 = minSR(Q1, Q2, Q3);
1223
1224 if (mn1 > mx2) return 0;
1225 if (mn2 > mx1) return 0;
1226 return 1;
1227 }
1228 //=======================================================================
1229 //function : TriContact
1230 //purpose : This function checks if two triangles are in
1231 // contact or not, return 1 if yes, return 0
1232 // if no.
1233 //=======================================================================
TriContact(const IntPolyh_Point & P1,const IntPolyh_Point & P2,const IntPolyh_Point & P3,const IntPolyh_Point & Q1,const IntPolyh_Point & Q2,const IntPolyh_Point & Q3,Standard_Real & Angle) const1234 Standard_Integer IntPolyh_MaillageAffinage::TriContact
1235 (const IntPolyh_Point &P1,
1236 const IntPolyh_Point &P2,
1237 const IntPolyh_Point &P3,
1238 const IntPolyh_Point &Q1,
1239 const IntPolyh_Point &Q2,
1240 const IntPolyh_Point &Q3,
1241 Standard_Real &Angle) const
1242 {
1243 /**
1244 The first triangle is (p1,p2,p3). The other is (q1,q2,q3).
1245 The edges are (e1,e2,e3) and (f1,f2,f3).
1246 The normals are n1 and m1
1247 Outwards are (g1,g2,g3) and (h1,h2,h3).*/
1248
1249 if(maxSR(P1.X(),P2.X(),P3.X())<minSR(Q1.X(),Q2.X(),Q3.X())) return(0);
1250 if(maxSR(P1.Y(),P2.Y(),P3.Y())<minSR(Q1.Y(),Q2.Y(),Q3.Y())) return(0);
1251 if(maxSR(P1.Z(),P2.Z(),P3.Z())<minSR(Q1.Z(),Q2.Z(),Q3.Z())) return(0);
1252
1253 if(minSR(P1.X(),P2.X(),P3.X())>maxSR(Q1.X(),Q2.X(),Q3.X())) return(0);
1254 if(minSR(P1.Y(),P2.Y(),P3.Y())>maxSR(Q1.Y(),Q2.Y(),Q3.Y())) return(0);
1255 if(minSR(P1.Z(),P2.Z(),P3.Z())>maxSR(Q1.Z(),Q2.Z(),Q3.Z())) return(0);
1256
1257 IntPolyh_Point p1, p2, p3;
1258 IntPolyh_Point q1, q2, q3;
1259 IntPolyh_Point e1, e2, e3;
1260 IntPolyh_Point f1, f2, f3;
1261 IntPolyh_Point g1, g2, g3;
1262 IntPolyh_Point h1, h2, h3;
1263 IntPolyh_Point n1, m1;
1264 IntPolyh_Point z;
1265
1266 IntPolyh_Point ef11, ef12, ef13;
1267 IntPolyh_Point ef21, ef22, ef23;
1268 IntPolyh_Point ef31, ef32, ef33;
1269
1270 z.SetX(0.0); z.SetY(0.0); z.SetZ(0.0);
1271
1272
1273 p1.SetX(P1.X() - P1.X()); p1.SetY(P1.Y() - P1.Y()); p1.SetZ(P1.Z() - P1.Z());
1274 p2.SetX(P2.X() - P1.X()); p2.SetY(P2.Y() - P1.Y()); p2.SetZ(P2.Z() - P1.Z());
1275 p3.SetX(P3.X() - P1.X()); p3.SetY(P3.Y() - P1.Y()); p3.SetZ(P3.Z() - P1.Z());
1276
1277 q1.SetX(Q1.X() - P1.X()); q1.SetY(Q1.Y() - P1.Y()); q1.SetZ(Q1.Z() - P1.Z());
1278 q2.SetX(Q2.X() - P1.X()); q2.SetY(Q2.Y() - P1.Y()); q2.SetZ(Q2.Z() - P1.Z());
1279 q3.SetX(Q3.X() - P1.X()); q3.SetY(Q3.Y() - P1.Y()); q3.SetZ(Q3.Z() - P1.Z());
1280
1281 e1.SetX(p2.X() - p1.X()); e1.SetY(p2.Y() - p1.Y()); e1.SetZ(p2.Z() - p1.Z());
1282 e2.SetX(p3.X() - p2.X()); e2.SetY(p3.Y() - p2.Y()); e2.SetZ(p3.Z() - p2.Z());
1283 e3.SetX(p1.X() - p3.X()); e3.SetY(p1.Y() - p3.Y()); e3.SetZ(p1.Z() - p3.Z());
1284
1285 f1.SetX(q2.X() - q1.X()); f1.SetY(q2.Y() - q1.Y()); f1.SetZ(q2.Z() - q1.Z());
1286 f2.SetX(q3.X() - q2.X()); f2.SetY(q3.Y() - q2.Y()); f2.SetZ(q3.Z() - q2.Z());
1287 f3.SetX(q1.X() - q3.X()); f3.SetY(q1.Y() - q3.Y()); f3.SetZ(q1.Z() - q3.Z());
1288
1289 n1.Cross(e1, e2); //normal to the first triangle
1290 m1.Cross(f1, f2); //normal to the second triangle
1291
1292 g1.Cross(e1, n1);
1293 g2.Cross(e2, n1);
1294 g3.Cross(e3, n1);
1295 h1.Cross(f1, m1);
1296 h2.Cross(f2, m1);
1297 h3.Cross(f3, m1);
1298
1299 ef11.Cross(e1, f1);
1300 ef12.Cross(e1, f2);
1301 ef13.Cross(e1, f3);
1302 ef21.Cross(e2, f1);
1303 ef22.Cross(e2, f2);
1304 ef23.Cross(e2, f3);
1305 ef31.Cross(e3, f1);
1306 ef32.Cross(e3, f2);
1307 ef33.Cross(e3, f3);
1308
1309 // Now the testing is done
1310
1311 if (!project6(n1, p1, p2, p3, q1, q2, q3)) return 0; //T2 is not higher or lower than T1
1312 if (!project6(m1, p1, p2, p3, q1, q2, q3)) return 0; //T1 is not higher of lower than T2
1313
1314 if (!project6(ef11, p1, p2, p3, q1, q2, q3)) return 0;
1315 if (!project6(ef12, p1, p2, p3, q1, q2, q3)) return 0;
1316 if (!project6(ef13, p1, p2, p3, q1, q2, q3)) return 0;
1317 if (!project6(ef21, p1, p2, p3, q1, q2, q3)) return 0;
1318 if (!project6(ef22, p1, p2, p3, q1, q2, q3)) return 0;
1319 if (!project6(ef23, p1, p2, p3, q1, q2, q3)) return 0;
1320 if (!project6(ef31, p1, p2, p3, q1, q2, q3)) return 0;
1321 if (!project6(ef32, p1, p2, p3, q1, q2, q3)) return 0;
1322 if (!project6(ef33, p1, p2, p3, q1, q2, q3)) return 0;
1323
1324 if (!project6(g1, p1, p2, p3, q1, q2, q3)) return 0; //T2 is outside of T1 in the plane of T1
1325 if (!project6(g2, p1, p2, p3, q1, q2, q3)) return 0; //T2 is outside of T1 in the plane of T1
1326 if (!project6(g3, p1, p2, p3, q1, q2, q3)) return 0; //T2 is outside of T1 in the plane of T1
1327 if (!project6(h1, p1, p2, p3, q1, q2, q3)) return 0; //T1 is outside of T2 in the plane of T2
1328 if (!project6(h2, p1, p2, p3, q1, q2, q3)) return 0; //T1 is outside of T2 in the plane of T2
1329 if (!project6(h3, p1, p2, p3, q1, q2, q3)) return 0; //T1 is outside of T2 in the plane of T2
1330
1331 //Calculation of cosinus angle between two normals
1332 Standard_Real SqModn1=-1.0;
1333 Standard_Real SqModm1=-1.0;
1334 SqModn1=n1.SquareModulus();
1335 if (SqModn1>SquareMyConfusionPrecision){
1336 SqModm1=m1.SquareModulus();
1337 }
1338 if (SqModm1>SquareMyConfusionPrecision) {
1339 Angle=(n1.Dot(m1))/(sqrt(SqModn1)*sqrt(SqModm1));
1340 }
1341 return 1;
1342 }
1343 //=======================================================================
1344 //function : TestNbPoints
1345 //purpose : This function is used by StartingPointsResearch() to control
1346 // the number of points found keep the result in conformity (1 or 2 points)
1347 // void TestNbPoints(const Standard_Integer TriSurfID,
1348 //=======================================================================
TestNbPoints(const Standard_Integer,Standard_Integer & NbPoints,Standard_Integer & NbPointsTotal,const IntPolyh_StartPoint & Pt1,const IntPolyh_StartPoint & Pt2,IntPolyh_StartPoint & SP1,IntPolyh_StartPoint & SP2)1349 void TestNbPoints(const Standard_Integer ,
1350 Standard_Integer &NbPoints,
1351 Standard_Integer &NbPointsTotal,
1352 const IntPolyh_StartPoint &Pt1,
1353 const IntPolyh_StartPoint &Pt2,
1354 IntPolyh_StartPoint &SP1,
1355 IntPolyh_StartPoint &SP2)
1356 {
1357 // already checked in TriangleEdgeContact
1358 // if( (NbPoints==2)&&(Pt1.CheckSameSP(Pt2)) ) NbPoints=1;
1359
1360 if(NbPoints>2) {
1361
1362 }
1363 else {
1364 if ( (NbPoints==1)&&(NbPointsTotal==0) ) {
1365 SP1=Pt1;
1366 NbPointsTotal=1;
1367 }
1368 else if ( (NbPoints==1)&&(NbPointsTotal==1) ) {
1369 if(Pt1.CheckSameSP(SP1)!=1) {
1370 SP2=Pt1;
1371 NbPointsTotal=2;
1372 }
1373 }
1374 else if( (NbPoints==1)&&(NbPointsTotal==2) ) {
1375 if ( (SP1.CheckSameSP(Pt1))||(SP2.CheckSameSP(Pt1)) )
1376 NbPointsTotal=2;
1377 else NbPointsTotal=3;
1378 }
1379 else if( (NbPoints==2)&&(NbPointsTotal==0) ) {
1380 SP1=Pt1;
1381 SP2=Pt2;
1382 NbPointsTotal=2;
1383 }
1384 else if( (NbPoints==2)&&(NbPointsTotal==1) ) {//there is also Pt1 != Pt2
1385 if(SP1.CheckSameSP(Pt1)) {
1386 SP2=Pt2;
1387 NbPointsTotal=2;
1388 }
1389 else if (SP1.CheckSameSP(Pt2)) {
1390 SP2=Pt1;
1391 NbPointsTotal=2;
1392 }
1393 else NbPointsTotal=3;///case SP1!=Pt1 && SP1!=Pt2!
1394 }
1395 else if( (NbPoints==2)&&(NbPointsTotal==2) ) {//there is also SP1!=SP2
1396 if( (SP1.CheckSameSP(Pt1))||(SP1.CheckSameSP(Pt2)) ) {
1397 if( (SP2.CheckSameSP(Pt1))||(SP2.CheckSameSP(Pt2)) )
1398 NbPointsTotal=2;
1399 else NbPointsTotal=3;
1400 }
1401 else NbPointsTotal=3;
1402 }
1403 }
1404 }
1405 //=======================================================================
1406 //function : StartingPointsResearch
1407 //purpose : From two triangles compute intersection points.
1408 // If I found more than two intersection points
1409 // it means that those triangle are coplanar
1410 //=======================================================================
StartingPointsResearch(const Standard_Integer T1,const Standard_Integer T2,IntPolyh_StartPoint & SP1,IntPolyh_StartPoint & SP2) const1411 Standard_Integer IntPolyh_MaillageAffinage::StartingPointsResearch
1412 (const Standard_Integer T1,
1413 const Standard_Integer T2,
1414 IntPolyh_StartPoint &SP1,
1415 IntPolyh_StartPoint &SP2) const
1416 {
1417 const IntPolyh_Triangle &Tri1=TTriangles1[T1];
1418 const IntPolyh_Triangle &Tri2=TTriangles2[T2];
1419
1420 const IntPolyh_Point &P1=TPoints1[Tri1.FirstPoint()];
1421 const IntPolyh_Point &P2=TPoints1[Tri1.SecondPoint()];
1422 const IntPolyh_Point &P3=TPoints1[Tri1.ThirdPoint()];
1423 const IntPolyh_Point &Q1=TPoints2[Tri2.FirstPoint()];
1424 const IntPolyh_Point &Q2=TPoints2[Tri2.SecondPoint()];
1425 const IntPolyh_Point &Q3=TPoints2[Tri2.ThirdPoint()];
1426
1427
1428
1429 /* The first triangle is (p1,p2,p3). The other is (q1,q2,q3).
1430 The sides are (e1,e2,e3) and (f1,f2,f3).
1431 The normals are n1 and m1*/
1432
1433 const IntPolyh_Point e1=P2-P1;
1434 const IntPolyh_Point e2=P3-P2;
1435 const IntPolyh_Point e3=P1-P3;
1436
1437 const IntPolyh_Point f1=Q2-Q1;
1438 const IntPolyh_Point f2=Q3-Q2;
1439 const IntPolyh_Point f3=Q1-Q3;
1440
1441
1442 IntPolyh_Point nn1,mm1;
1443 nn1.Cross(e1, e2); //normal to the first triangle
1444 mm1.Cross(f1, f2); //normal to the second triangle
1445
1446 Standard_Real nn1modulus, mm1modulus;
1447 nn1modulus=sqrt(nn1.SquareModulus());
1448 mm1modulus=sqrt(mm1.SquareModulus());
1449
1450 //-------------------------------------------------
1451 ///calculation of intersections points between triangles
1452 //-------------------------------------------------
1453 Standard_Integer NbPoints=0;
1454 Standard_Integer NbPointsTotal=0;
1455
1456
1457 ///check T1 normal
1458 if(Abs(nn1modulus)<MyConfusionPrecision) {//10.0e-20){
1459
1460 }
1461 else {
1462 const IntPolyh_Point n1=nn1.Divide(nn1modulus);
1463 ///T2 edges with T1
1464 if(NbPointsTotal<3) {
1465 IntPolyh_StartPoint Pt1,Pt2;
1466 NbPoints=TriangleEdgeContact(1,1,Tri1,Tri2,P1,P2,P3,e1,e2,e3,Q1,Q2,f1,n1,Pt1,Pt2);
1467 TestNbPoints(1,NbPoints,NbPointsTotal,Pt1,Pt2,SP1,SP2);
1468 }
1469
1470 if(NbPointsTotal<3) {
1471 IntPolyh_StartPoint Pt1,Pt2;
1472 NbPoints=TriangleEdgeContact(1,2,Tri1,Tri2,P1,P2,P3,e1,e2,e3,Q2,Q3,f2,n1,Pt1,Pt2);
1473 TestNbPoints(1,NbPoints,NbPointsTotal,Pt1,Pt2,SP1,SP2);
1474 }
1475
1476 if(NbPointsTotal<3) {
1477 IntPolyh_StartPoint Pt1,Pt2;
1478 NbPoints=TriangleEdgeContact(1,3,Tri1,Tri2,P1,P2,P3,e1,e2,e3,Q3,Q1,f3,n1,Pt1,Pt2);
1479 TestNbPoints(1,NbPoints,NbPointsTotal,Pt1,Pt2,SP1,SP2);
1480 }
1481 }
1482
1483 ///check T2 normal
1484 if(Abs(mm1modulus)<MyConfusionPrecision) {//10.0e-20){
1485
1486 }
1487 else {
1488 const IntPolyh_Point m1=mm1.Divide(mm1modulus);
1489 ///T1 edges with T2
1490 if(NbPointsTotal<3) {
1491 IntPolyh_StartPoint Pt1,Pt2;
1492 NbPoints=TriangleEdgeContact(2,1,Tri1,Tri2,Q1,Q2,Q3,f1,f2,f3,P1,P2,e1,m1,Pt1,Pt2);
1493 TestNbPoints(2,NbPoints,NbPointsTotal,Pt1,Pt2,SP1,SP2);
1494 }
1495
1496 if(NbPointsTotal<3) {
1497 IntPolyh_StartPoint Pt1,Pt2;
1498 NbPoints=TriangleEdgeContact(2,2,Tri1,Tri2,Q1,Q2,Q3,f1,f2,f3,P2,P3,e2,m1,Pt1,Pt2);
1499 TestNbPoints(2,NbPoints,NbPointsTotal,Pt1,Pt2,SP1,SP2);
1500 }
1501
1502 if(NbPointsTotal<3) {
1503 IntPolyh_StartPoint Pt1,Pt2;
1504 NbPoints=TriangleEdgeContact(2,3,Tri1,Tri2,Q1,Q2,Q3,f1,f2,f3,P3,P1,e3,m1,Pt1,Pt2);
1505 TestNbPoints(2,NbPoints,NbPointsTotal,Pt1,Pt2,SP1,SP2);
1506 }
1507 }
1508 // no need already checked in TestNbPoints
1509 // if( (NbPointsTotal==2)&&(SP1.CheckSameSP(SP2)) ) {
1510 // NbPointsTotal=1;
1511 //SP1.SetCoupleValue(T1,T2);
1512 // }
1513 // else
1514 if(NbPointsTotal==2) {
1515 SP1.SetCoupleValue(T1,T2);
1516 SP2.SetCoupleValue(T1,T2);
1517 }
1518 else if(NbPointsTotal==1)
1519 SP1.SetCoupleValue(T1,T2);
1520 else if(NbPointsTotal==3)
1521 SP1.SetCoupleValue(T1,T2);
1522
1523 return (NbPointsTotal);
1524 }
1525 //=======================================================================
1526 //function : NextStartingPointsResearch
1527 //purpose : from two triangles and an intersection point I
1528 // search the other point (if it exist).
1529 // This function is used by StartPointChain
1530 //=======================================================================
NextStartingPointsResearch(const Standard_Integer T1,const Standard_Integer T2,const IntPolyh_StartPoint & SPInit,IntPolyh_StartPoint & SPNext) const1531 Standard_Integer IntPolyh_MaillageAffinage::NextStartingPointsResearch
1532 (const Standard_Integer T1,
1533 const Standard_Integer T2,
1534 const IntPolyh_StartPoint &SPInit,
1535 IntPolyh_StartPoint &SPNext) const
1536 {
1537 Standard_Integer NbPointsTotal=0;
1538 Standard_Integer EdgeInit1=SPInit.E1();
1539 Standard_Integer EdgeInit2=SPInit.E2();
1540 if( (T1<0)||(T2<0) ) NbPointsTotal=0;
1541 else {
1542
1543 const IntPolyh_Triangle &Tri1=TTriangles1[T1];
1544 const IntPolyh_Triangle &Tri2=TTriangles2[T2];
1545
1546 const IntPolyh_Point &P1=TPoints1[Tri1.FirstPoint()];
1547 const IntPolyh_Point &P2=TPoints1[Tri1.SecondPoint()];
1548 const IntPolyh_Point &P3=TPoints1[Tri1.ThirdPoint()];
1549 const IntPolyh_Point &Q1=TPoints2[Tri2.FirstPoint()];
1550 const IntPolyh_Point &Q2=TPoints2[Tri2.SecondPoint()];
1551 const IntPolyh_Point &Q3=TPoints2[Tri2.ThirdPoint()];
1552
1553 /* The first triangle is (p1,p2,p3). The other is (q1,q2,q3).
1554 The edges are (e1,e2,e3) and (f1,f2,f3).
1555 The normals are n1 and m1*/
1556
1557 const IntPolyh_Point e1=P2-P1;
1558 const IntPolyh_Point e2=P3-P2;
1559 const IntPolyh_Point e3=P1-P3;
1560
1561 const IntPolyh_Point f1=Q2-Q1;
1562 const IntPolyh_Point f2=Q3-Q2;
1563 const IntPolyh_Point f3=Q1-Q3;
1564
1565 IntPolyh_Point nn1,mm1;
1566 nn1.Cross(e1, e2); //normal to the first triangle
1567 mm1.Cross(f1, f2); //normal to the second triangle
1568
1569 Standard_Real nn1modulus, mm1modulus;
1570 nn1modulus=sqrt(nn1.SquareModulus());
1571 mm1modulus=sqrt(mm1.SquareModulus());
1572
1573 //-------------------------------------------------
1574 ///calculation of intersections points between triangles
1575 //-------------------------------------------------
1576
1577 Standard_Integer NbPoints=0;
1578 IntPolyh_StartPoint SP1,SP2;
1579
1580 ///check T1 normal
1581 if(Abs(nn1modulus)<MyConfusionPrecision) {//10.0e-20){
1582
1583 }
1584 else {
1585 const IntPolyh_Point n1=nn1.Divide(nn1modulus);
1586 ///T2 edges with T1
1587 if( (NbPointsTotal<3)&&(EdgeInit2!=Tri2.FirstEdge()) ) {
1588 IntPolyh_StartPoint Pt1,Pt2;
1589 NbPoints=TriangleEdgeContact(1,1,Tri1,Tri2,P1,P2,P3,e1,e2,e3,Q1,Q2,f1,n1,Pt1,Pt2);
1590 TestNbPoints(1,NbPoints,NbPointsTotal,Pt1,Pt2,SP1,SP2);
1591 }
1592
1593 if( (NbPointsTotal<3)&&(EdgeInit2!=Tri2.SecondEdge()) ) {
1594 IntPolyh_StartPoint Pt1,Pt2;
1595 NbPoints=TriangleEdgeContact(1,2,Tri1,Tri2,P1,P2,P3,e1,e2,e3,Q2,Q3,f2,n1,Pt1,Pt2);
1596 TestNbPoints(1,NbPoints,NbPointsTotal,Pt1,Pt2,SP1,SP2);
1597 }
1598
1599 if( (NbPointsTotal<3)&&(EdgeInit2!=Tri2.ThirdEdge()) ) {
1600 IntPolyh_StartPoint Pt1,Pt2;
1601 NbPoints=TriangleEdgeContact(1,3,Tri1,Tri2,P1,P2,P3,e1,e2,e3,Q3,Q1,f3,n1,Pt1,Pt2);
1602 TestNbPoints(1,NbPoints,NbPointsTotal,Pt1,Pt2,SP1,SP2);
1603 }
1604 }
1605 ///check T2 normal
1606 if(Abs(mm1modulus)<MyConfusionPrecision) {//10.0e-20){
1607
1608 }
1609 else {
1610 const IntPolyh_Point m1=mm1.Divide(mm1modulus);
1611 ///T1 edges with T2
1612 if( (NbPointsTotal<3)&&(EdgeInit1!=Tri1.FirstEdge()) ) {
1613 IntPolyh_StartPoint Pt1,Pt2;
1614 NbPoints=TriangleEdgeContact(2,1,Tri1,Tri2,Q1,Q2,Q3,f1,f2,f3,P1,P2,e1,m1,Pt1,Pt2);
1615 TestNbPoints(2,NbPoints,NbPointsTotal,Pt1,Pt2,SP1,SP2);
1616 }
1617
1618 if( (NbPointsTotal<3)&&(EdgeInit1!=Tri1.SecondEdge()) ) {
1619 IntPolyh_StartPoint Pt1,Pt2;
1620 NbPoints=TriangleEdgeContact(2,2,Tri1,Tri2,Q1,Q2,Q3,f1,f2,f3,P2,P3,e2,m1,Pt1,Pt2);
1621 TestNbPoints(2,NbPoints,NbPointsTotal,Pt1,Pt2,SP1,SP2);
1622 }
1623
1624 if( (NbPointsTotal<3)&&(EdgeInit1!=Tri1.ThirdEdge()) ) {
1625 IntPolyh_StartPoint Pt1,Pt2;
1626 NbPoints=TriangleEdgeContact(2,3,Tri1,Tri2,Q1,Q2,Q3,f1,f2,f3,P3,P1,e3,m1,Pt1,Pt2);
1627 TestNbPoints(2,NbPoints,NbPointsTotal,Pt1,Pt2,SP1,SP2);
1628 }
1629 }
1630
1631 if (NbPointsTotal==1) {
1632 if(SP1.CheckSameSP(SPInit))
1633 NbPointsTotal=0;
1634 else {
1635 SPNext=SP1;
1636 }
1637 }
1638 else if( (NbPointsTotal==2)&&(SP1.CheckSameSP(SPInit)) ) {
1639 NbPointsTotal=1;//SP1 et SPInit sont identiques
1640 SPNext=SP2;
1641 }
1642 else if( (NbPointsTotal==2)&&(SP2.CheckSameSP(SPInit)) ) {
1643 NbPointsTotal=1;//SP2 et SPInit sont identiques
1644 SPNext=SP1;
1645 }
1646
1647 else if(NbPointsTotal>1) {
1648
1649 }
1650 }
1651 SPNext.SetCoupleValue(T1,T2);
1652 return (NbPointsTotal);
1653 }
1654 //=======================================================================
1655 //function : CalculPtsInterTriEdgeCoplanaires
1656 //purpose :
1657 //=======================================================================
CalculPtsInterTriEdgeCoplanaires(const Standard_Integer TriSurfID,const IntPolyh_Point & NormaleTri,const IntPolyh_Triangle & Tri1,const IntPolyh_Triangle & Tri2,const IntPolyh_Point & PE1,const IntPolyh_Point & PE2,const IntPolyh_Point & Edge,const Standard_Integer EdgeIndex,const IntPolyh_Point & PT1,const IntPolyh_Point & PT2,const IntPolyh_Point & Cote,const Standard_Integer CoteIndex,IntPolyh_StartPoint & SP1,IntPolyh_StartPoint & SP2,Standard_Integer & NbPoints)1658 void CalculPtsInterTriEdgeCoplanaires(const Standard_Integer TriSurfID,
1659 const IntPolyh_Point &NormaleTri,
1660 const IntPolyh_Triangle &Tri1,
1661 const IntPolyh_Triangle &Tri2,
1662 const IntPolyh_Point &PE1,
1663 const IntPolyh_Point &PE2,
1664 const IntPolyh_Point &Edge,
1665 const Standard_Integer EdgeIndex,
1666 const IntPolyh_Point &PT1,
1667 const IntPolyh_Point &PT2,
1668 const IntPolyh_Point &Cote,
1669 const Standard_Integer CoteIndex,
1670 IntPolyh_StartPoint &SP1,
1671 IntPolyh_StartPoint &SP2,
1672 Standard_Integer &NbPoints)
1673 {
1674 Standard_Real aDE, aDC;
1675 //
1676 gp_Vec aVE(Edge.X(), Edge.Y(), Edge.Z());
1677 gp_Vec aVC(Cote.X(), Cote.Y(), Cote.Z());
1678 //
1679 aDE = aVE.SquareMagnitude();
1680 aDC = aVC.SquareMagnitude();
1681 //
1682 if (aDE > SquareMyConfusionPrecision) {
1683 aVE.Divide(aDE);
1684 }
1685 //
1686 if (aDC > SquareMyConfusionPrecision) {
1687 aVC.Divide(aDC);
1688 }
1689 //
1690 if (!aVE.IsParallel(aVC, MyConfusionPrecision)) {
1691 ///Edge and side are not parallel
1692 IntPolyh_Point Per;
1693 Per.Cross(NormaleTri,Cote);
1694 Standard_Real p1p = Per.Dot(PE1);
1695 Standard_Real p2p = Per.Dot(PE2);
1696 Standard_Real p0p = Per.Dot(PT1);
1697 ///The edge are PT1 are projected on the perpendicular of the side in the plane of the triangle
1698 if ( ( ((p1p>=p0p)&&(p0p>=p2p) )||( (p1p<=p0p)&&(p0p<=p2p) )) && (Abs(p1p-p2p) > MyConfusionPrecision)) {
1699 Standard_Real lambda=(p1p-p0p)/(p1p-p2p);
1700 if (lambda<-MyConfusionPrecision) {
1701
1702 }
1703 IntPolyh_Point PIE;
1704 if (Abs(lambda)<MyConfusionPrecision)//lambda=0
1705 PIE=PE1;
1706 else if (Abs(lambda)>1.0-MyConfusionPrecision)//lambda=1
1707 PIE=PE2;
1708 else
1709 PIE=PE1+Edge*lambda;
1710
1711 Standard_Real alpha=RealLast();
1712 if(Cote.X()!=0) alpha=(PIE.X()-PT1.X())/Cote.X();
1713 else if (Cote.Y()!=0) alpha=(PIE.Y()-PT1.Y())/Cote.Y();
1714 else if (Cote.Z()!=0) alpha=(PIE.Z()-PT1.Z())/Cote.Z();
1715 else {
1716
1717 }
1718
1719 if (alpha<-MyConfusionPrecision) {
1720
1721 }
1722 else {
1723 if (NbPoints==0) {
1724 SP1.SetXYZ(PIE.X(),PIE.Y(),PIE.Z());
1725 if (TriSurfID==1) {
1726 if(Abs(alpha)<MyConfusionPrecision) {//alpha=0
1727 SP1.SetUV1(PT1.U(),PT1.V());
1728 SP1.SetUV1(PIE.U(),PIE.V());
1729 SP1.SetEdge1(-1);
1730 }
1731 if(Abs(alpha)>1.0-MyConfusionPrecision) {//alpha=1
1732 SP1.SetUV1(PT2.U(),PT2.V());
1733 SP1.SetUV1(PIE.U(),PIE.V());
1734 SP1.SetEdge1(-1);
1735 }
1736 else {
1737 SP1.SetUV1(PT1.U()+Cote.U()*alpha,PT1.V()+Cote.V()*alpha);
1738 SP1.SetUV2(PIE.U(),PIE.V());
1739 SP1.SetEdge1(Tri1.GetEdgeNumber(CoteIndex));
1740 if (Tri1.GetEdgeOrientation(CoteIndex)>0) SP1.SetLambda1(alpha);
1741 else SP1.SetLambda1(1.0-alpha);
1742 }
1743 NbPoints++;
1744 }
1745 else if (TriSurfID==2) {
1746 if(Abs(alpha)<MyConfusionPrecision) {//alpha=0
1747 SP1.SetUV1(PT1.U(),PT1.V());
1748 SP1.SetUV1(PIE.U(),PIE.V());
1749 SP1.SetEdge2(-1);
1750 }
1751 if(Abs(alpha)>1.0-MyConfusionPrecision) {//alpha=1
1752 SP1.SetUV1(PT2.U(),PT2.V());
1753 SP1.SetUV1(PIE.U(),PIE.V());
1754 SP1.SetEdge2(-1);
1755 }
1756 else {
1757 SP1.SetUV1(PIE.U(),PIE.V());
1758 SP1.SetUV2(PT1.U()+Cote.U()*alpha,PT1.V()+Cote.V()*alpha);
1759 SP1.SetEdge2(Tri2.GetEdgeNumber(CoteIndex));
1760 if (Tri2.GetEdgeOrientation(CoteIndex)>0) SP1.SetLambda2(alpha);
1761 else SP1.SetLambda2(1.0-alpha);
1762 }
1763 NbPoints++;
1764 }
1765 else {
1766
1767 }
1768 }
1769
1770 else if (NbPoints==1) {
1771 SP2.SetXYZ(PIE.X(),PIE.Y(),PIE.Z());
1772 if (TriSurfID==1) {
1773 if(Abs(alpha)<MyConfusionPrecision) {//alpha=0
1774 SP2.SetUV1(PT1.U(),PT1.V());
1775 SP2.SetUV1(PIE.U(),PIE.V());
1776 SP2.SetEdge1(-1);
1777 }
1778 if(Abs(alpha)>1.0-MyConfusionPrecision) {//alpha=1
1779 SP2.SetUV1(PT2.U(),PT2.V());
1780 SP2.SetUV1(PIE.U(),PIE.V());
1781 SP2.SetEdge1(-1);
1782 }
1783 else {
1784 SP2.SetUV1(PT1.U()+Cote.U()*alpha,PT1.V()+Cote.V()*alpha);
1785 SP2.SetUV2(PIE.U(),PIE.V());
1786 SP2.SetEdge1(Tri1.GetEdgeNumber(CoteIndex));
1787 if (Tri1.GetEdgeOrientation(CoteIndex)>0) SP2.SetLambda1(alpha);
1788 else SP2.SetLambda1(1.0-alpha);
1789 }
1790 NbPoints++;
1791 }
1792 else if (TriSurfID==2) {
1793 if(Abs(alpha)<MyConfusionPrecision) {//alpha=0
1794 SP2.SetUV1(PT1.U(),PT1.V());
1795 SP2.SetUV1(PIE.U(),PIE.V());
1796 SP2.SetEdge2(-1);
1797 }
1798 if(Abs(alpha)>1.0-MyConfusionPrecision) {//alpha=1
1799 SP2.SetUV1(PT2.U(),PT2.V());
1800 SP2.SetUV1(PIE.U(),PIE.V());
1801 SP2.SetEdge2(-1);
1802 }
1803 else {
1804 SP2.SetUV1(PIE.U(),PIE.V());
1805 SP2.SetUV2(PT1.U()+Cote.U()*alpha,PT1.V()+Cote.V()*alpha);
1806 SP2.SetEdge2(Tri2.GetEdgeNumber(CoteIndex));
1807 if (Tri1.GetEdgeOrientation(CoteIndex)>0) SP2.SetLambda2(alpha);
1808 else SP2.SetLambda2(1.0-alpha);
1809 }
1810 NbPoints++;
1811 }
1812 else {
1813
1814 }
1815 }
1816
1817 else if( (NbPoints>2)||(NbPoints<0) ) {
1818
1819 }
1820 }
1821 }
1822 }
1823 else {
1824 //Side and Edge are parallel, with previous
1825 //rejections they are at the same side
1826 //The points are projected on that side
1827 Standard_Real pe1p= Cote.Dot(PE1);
1828 Standard_Real pe2p= Cote.Dot(PE2);
1829 Standard_Real pt1p= Cote.Dot(PT1);
1830 Standard_Real pt2p= Cote.Dot(PT2);
1831 Standard_Real lambda1=0., lambda2=0., alpha1=0., alpha2=0.;
1832 IntPolyh_Point PEP1,PTP1,PEP2,PTP2;
1833
1834 if (pe1p>pe2p) {
1835 if ( (pt1p<pe1p)&&(pe1p<=pt2p) ) {
1836 lambda1=0.0;
1837 PEP1=PE1;
1838 alpha1=((pe1p-pt1p)/(pt2p-pt1p));
1839 PTP1=PT1+Cote*alpha1;
1840 NbPoints=1;
1841 if (pt1p<=pe2p) {
1842 lambda2=1.0;
1843 PEP2=PE2;
1844 alpha2=((pe2p-pt1p)/(pt2p-pt1p));
1845 PTP2=PT1+Cote*alpha2;
1846 NbPoints=2;
1847 }
1848 else {
1849 lambda2=((pt1p-pe1p)/(pe2p-pe1p));
1850 PEP2=PE1+Edge*lambda2;
1851 alpha2=0.0;
1852 PTP2=PT1;
1853 NbPoints=2;
1854 }
1855 }
1856 else if( (pt2p<pe1p)&&(pe1p<=pt1p) ) {
1857 lambda1=0.0;
1858 PEP1=PE1;
1859 alpha1=((pt1p-pe1p)/(pt1p-pt2p));
1860 PTP1=PT1+Cote*alpha1;
1861 NbPoints=1;
1862 if (pt2p<=pe2p) {
1863 lambda2=1.0;
1864 PEP2=PE2;
1865 alpha2=((pe2p-pt1p)/(pt2p-pt1p));
1866 PTP2=PT1+Cote*alpha2;
1867 NbPoints=2;
1868 }
1869 else {
1870 lambda2=((pt2p-pe1p)/(pe2p-pe1p));
1871 PEP2=PE1+Edge*lambda2;
1872 alpha2=1.0;
1873 PTP2=PT2;
1874 NbPoints=2;
1875 }
1876 }
1877 }
1878
1879 if (pe1p<pe2p) {
1880 if ( (pt1p<pe2p)&&(pe2p<=pt2p) ) {
1881 lambda1=1.0;
1882 PEP1=PE2;
1883 alpha1=((pe2p-pt1p)/(pt2p-pt1p));
1884 PTP1=PT1+Cote*alpha1;
1885 NbPoints=1;
1886 if (pt1p<=pe1p) {
1887 lambda2=0.0;
1888 PEP2=PE1;
1889 alpha2=((pe1p-pt1p)/(pt2p-pt1p));
1890 PTP2=PT1+Cote*alpha2;
1891 NbPoints=2;
1892 }
1893 else {
1894 lambda2=((pt1p-pe1p)/(pe2p-pe1p));
1895 PEP2=PE2+Edge*lambda2;
1896 alpha2=0.0;
1897 PTP2=PT1;
1898 NbPoints=2;
1899 }
1900 }
1901 else if( (pt2p<pe2p)&&(pe2p<=pt1p) ) {
1902 lambda1=1.0;
1903 PEP1=PE2;
1904 alpha1=((pt1p-pe2p)/(pt1p-pt2p));
1905 PTP1=PT1+Cote*alpha1;
1906 NbPoints=1;
1907 if (pt2p<=pe1p) {
1908 lambda2=0.0;
1909 PEP2=PE1;
1910 alpha2=((pe1p-pt1p)/(pt2p-pt1p));
1911 PTP2=PT1+Cote*alpha2;
1912 NbPoints=2;
1913 }
1914 else {
1915 lambda2=((pt2p-pe1p)/(pe2p-pe1p));
1916 PEP2=PE1+Edge*lambda2;
1917 alpha2=1.0;
1918 PTP2=PT2;
1919 NbPoints=2;
1920 }
1921 }
1922 }
1923
1924 if (NbPoints!=0) {
1925 SP1.SetXYZ(PEP1.X(),PEP1.Y(),PEP1.Z());
1926 if (TriSurfID==1) {///cote appartient a Tri1
1927 SP1.SetUV1(PTP1.U(),PTP1.V());
1928 SP1.SetUV2(PEP1.U(),PEP1.V());
1929 SP1.SetEdge1(Tri1.GetEdgeNumber(CoteIndex));
1930
1931 if(Tri1.GetEdgeOrientation(CoteIndex)>0) SP1.SetLambda1(alpha1);
1932 else SP1.SetLambda1(1.0-alpha1);
1933
1934 if(Tri2.GetEdgeOrientation(EdgeIndex)>0) SP1.SetLambda2(lambda1);
1935 else SP1.SetLambda2(1.0-lambda1);
1936 }
1937 if (TriSurfID==2) {///cote appartient a Tri2
1938 SP1.SetUV1(PEP1.U(),PTP1.V());
1939 SP1.SetUV2(PTP1.U(),PEP1.V());
1940 SP1.SetEdge2(Tri2.GetEdgeNumber(CoteIndex));
1941
1942 if(Tri2.GetEdgeOrientation(CoteIndex)>0) SP1.SetLambda1(alpha1);
1943 else SP1.SetLambda1(1.0-alpha1);
1944
1945 if(Tri1.GetEdgeOrientation(EdgeIndex)>0) SP1.SetLambda2(lambda1);
1946 else SP1.SetLambda2(1.0-lambda1);
1947 }
1948
1949 //It is checked if PEP1!=PEP2
1950 if ( (NbPoints==2)&&(Abs(PEP1.U()-PEP2.U())<MyConfusionPrecision)
1951 &&(Abs(PEP1.V()-PEP2.V())<MyConfusionPrecision) ) NbPoints=1;
1952 if (NbPoints==2) {
1953 SP2.SetXYZ(PEP2.X(),PEP2.Y(),PEP2.Z());
1954 if (TriSurfID==1) {
1955 SP2.SetUV1(PTP2.U(),PTP2.V());
1956 SP2.SetUV2(PEP2.U(),PEP2.V());
1957 SP2.SetEdge1(Tri1.GetEdgeNumber(CoteIndex));
1958
1959 if(Tri1.GetEdgeOrientation(CoteIndex)>0) SP2.SetLambda1(alpha1);
1960 else SP2.SetLambda1(1.0-alpha1);
1961
1962 if(Tri2.GetEdgeOrientation(EdgeIndex)>0) SP2.SetLambda2(lambda1);
1963 else SP2.SetLambda2(1.0-lambda1);
1964 }
1965 if (TriSurfID==2) {
1966 SP2.SetUV1(PEP2.U(),PTP2.V());
1967 SP2.SetUV2(PTP2.U(),PEP2.V());
1968 SP2.SetEdge2(Tri2.GetEdgeNumber(CoteIndex));
1969
1970 if(Tri1.GetEdgeOrientation(CoteIndex)>0) SP2.SetLambda1(alpha1);
1971 else SP2.SetLambda1(1.0-alpha1);
1972
1973 if(Tri2.GetEdgeOrientation(EdgeIndex)>0) SP2.SetLambda2(lambda1);
1974 else SP2.SetLambda2(1.0-lambda1);
1975 }
1976 }
1977 }
1978 }
1979 //Filter if the point is placed on top, the edge is set to -1
1980 if (NbPoints>0) {
1981 if(Abs(SP1.Lambda1())<MyConfusionPrecision)
1982 SP1.SetEdge1(-1);
1983 if(Abs(SP1.Lambda1()-1)<MyConfusionPrecision)
1984 SP1.SetEdge1(-1);
1985 if(Abs(SP1.Lambda2())<MyConfusionPrecision)
1986 SP1.SetEdge2(-1);
1987 if(Abs(SP1.Lambda2()-1)<MyConfusionPrecision)
1988 SP1.SetEdge2(-1);
1989 }
1990 if (NbPoints==2) {
1991 if(Abs(SP2.Lambda1())<MyConfusionPrecision)
1992 SP2.SetEdge1(-1);
1993 if(Abs(SP2.Lambda1()-1)<MyConfusionPrecision)
1994 SP2.SetEdge1(-1);
1995 if(Abs(SP2.Lambda2())<MyConfusionPrecision)
1996 SP2.SetEdge2(-1);
1997 if(Abs(SP2.Lambda2()-1)<MyConfusionPrecision)
1998 SP2.SetEdge2(-1);
1999 }
2000 }
2001
2002 //=======================================================================
2003 //function : TriangleEdgeContact
2004 //purpose :
2005 //=======================================================================
TriangleEdgeContact(const Standard_Integer TriSurfID,const Standard_Integer EdgeIndex,const IntPolyh_Triangle & Tri1,const IntPolyh_Triangle & Tri2,const IntPolyh_Point & PT1,const IntPolyh_Point & PT2,const IntPolyh_Point & PT3,const IntPolyh_Point & Cote12,const IntPolyh_Point & Cote23,const IntPolyh_Point & Cote31,const IntPolyh_Point & PE1,const IntPolyh_Point & PE2,const IntPolyh_Point & Edge,const IntPolyh_Point & NormaleT,IntPolyh_StartPoint & SP1,IntPolyh_StartPoint & SP2) const2006 Standard_Integer IntPolyh_MaillageAffinage::TriangleEdgeContact
2007 (const Standard_Integer TriSurfID,
2008 const Standard_Integer EdgeIndex,
2009 const IntPolyh_Triangle &Tri1,
2010 const IntPolyh_Triangle &Tri2,
2011 const IntPolyh_Point &PT1,
2012 const IntPolyh_Point &PT2,
2013 const IntPolyh_Point &PT3,
2014 const IntPolyh_Point &Cote12,
2015 const IntPolyh_Point &Cote23,
2016 const IntPolyh_Point &Cote31,
2017 const IntPolyh_Point &PE1,const IntPolyh_Point &PE2,
2018 const IntPolyh_Point &Edge,
2019 const IntPolyh_Point &NormaleT,
2020 IntPolyh_StartPoint &SP1,IntPolyh_StartPoint &SP2) const
2021 {
2022
2023 Standard_Real lambda =0., alpha =0., beta =0.;
2024
2025 //One of edges on which the point is located is known
2026 if (TriSurfID==1) {
2027 SP1.SetEdge2(Tri2.GetEdgeNumber(EdgeIndex));
2028 SP2.SetEdge2(Tri2.GetEdgeNumber(EdgeIndex));
2029 }
2030 else if (TriSurfID==2) {
2031 SP1.SetEdge1(Tri1.GetEdgeNumber(EdgeIndex));
2032 SP2.SetEdge1(Tri1.GetEdgeNumber(EdgeIndex));
2033 }
2034 else {
2035
2036 }
2037
2038 /**The edge is projected on the normal in the triangle if yes
2039 --> a free intersection (point I)--> Start point is found */
2040 Standard_Integer NbPoints=0;
2041 if(NormaleT.SquareModulus()==0) {
2042
2043 }
2044 else if( (Cote12.SquareModulus()==0)
2045 ||(Cote23.SquareModulus()==0)
2046 ||(Cote31.SquareModulus()==0) ) {
2047
2048 }
2049 else if(Edge.SquareModulus()==0) {
2050
2051 }
2052 else {
2053 Standard_Real pe1 = NormaleT.Dot(PE1);
2054 Standard_Real pe2 = NormaleT.Dot(PE2);
2055 Standard_Real pt1 = NormaleT.Dot(PT1);
2056
2057 // PE1I = lambda.Edge
2058
2059 if( (Abs(pe1-pt1)<MyConfusionPrecision)&&(Abs(pe2-pt1)<MyConfusionPrecision)) {
2060 //edge and triangle are coplanar (two contact points at maximum)
2061
2062
2063 //the tops of the triangle are projected on the perpendicular to the edge
2064 IntPolyh_Point PerpEdge;
2065 PerpEdge.Cross(NormaleT,Edge);
2066 Standard_Real pp1 = PerpEdge.Dot(PT1);
2067 Standard_Real pp2 = PerpEdge.Dot(PT2);
2068 Standard_Real pp3 = PerpEdge.Dot(PT3);
2069 Standard_Real ppe1 = PerpEdge.Dot(PE1);
2070
2071
2072 if ( (Abs(pp1-pp2)<MyConfusionPrecision)&&(Abs(pp1-pp3)<MyConfusionPrecision) ) {
2073
2074 }
2075 else {
2076 if ( ( (pp1>=ppe1)&&(pp2<=ppe1)&&(pp3<=ppe1) ) || ( (pp1<=ppe1)&&(pp2>=ppe1)&&(pp3>=ppe1) ) ){
2077 //there are two sides (common top PT1) that can cut the edge
2078
2079 //first side
2080 CalculPtsInterTriEdgeCoplanaires(TriSurfID,NormaleT,Tri1,Tri2,PE1,PE2,Edge,EdgeIndex,
2081 PT1,PT2,Cote12,1,SP1,SP2,NbPoints);
2082
2083 if ( (NbPoints>1)&&(Abs(SP2.U1()-SP1.U1())<MyConfusionPrecision)
2084 &&(Abs(SP1.V1()-SP2.V1())<MyConfusionPrecision) ) NbPoints=1;
2085
2086 //second side
2087 if (NbPoints<2) CalculPtsInterTriEdgeCoplanaires(TriSurfID,NormaleT,Tri1,Tri2,PE1,PE2,Edge,EdgeIndex,
2088 PT3,PT1,Cote31,3,SP1,SP2,NbPoints);
2089 }
2090
2091 if ( (NbPoints>1)&&(Abs(SP1.U1()-SP2.U1())<MyConfusionPrecision)
2092 &&(Abs(SP1.V2()-SP2.V1())<MyConfusionPrecision) ) NbPoints=1;
2093 if (NbPoints>=2) return(NbPoints);
2094
2095 else if ( ( ( (pp2>=ppe1)&&(pp1<=ppe1)&&(pp3<=ppe1) ) || ( (pp2<=ppe1)&&(pp1>=ppe1)&&(pp3>=ppe1) ) )
2096 && (NbPoints<2) ) {
2097 //there are two sides (common top PT2) that can cut the edge
2098
2099 //first side
2100 CalculPtsInterTriEdgeCoplanaires(TriSurfID,NormaleT,Tri1,Tri2,PE1,PE2,Edge,EdgeIndex,
2101 PT1,PT2,Cote12,1,SP1,SP2,NbPoints);
2102
2103 if ( (NbPoints>1)&&(Abs(SP2.U1()-SP1.U1())<MyConfusionPrecision)
2104 &&(Abs(SP1.V1()-SP2.V1())<MyConfusionPrecision) ) NbPoints=1;
2105
2106 //second side
2107 if(NbPoints<2) CalculPtsInterTriEdgeCoplanaires(TriSurfID,NormaleT,Tri1,Tri2,PE1,PE2,Edge,EdgeIndex,
2108 PT2,PT3,Cote23,2,SP1,SP2,NbPoints);
2109 }
2110 if ( (NbPoints>1)&&(Abs(SP2.U1()-SP1.U1())<MyConfusionPrecision)
2111 &&(Abs(SP1.V1()-SP2.V1())<MyConfusionPrecision) ) NbPoints=1;
2112 if (NbPoints>=2) return(NbPoints);
2113
2114 else if ( (( (pp3>=ppe1)&&(pp1<=ppe1)&&(pp2<=ppe1) ) || ( (pp3<=ppe1)&&(pp1>=ppe1)&&(pp2>=ppe1) ))
2115 && (NbPoints<2) ) {
2116 //there are two sides (common top PT3) that can cut the edge
2117
2118 //first side
2119 CalculPtsInterTriEdgeCoplanaires(TriSurfID,NormaleT,Tri1,Tri2,PE1,PE2,Edge,EdgeIndex,
2120 PT3,PT1,Cote31,3,SP1,SP2,NbPoints);
2121
2122 if ( (NbPoints>1)&&(Abs(SP2.U1()-SP1.U1())<MyConfusionPrecision)
2123 &&(Abs(SP1.V1()-SP2.V1())<MyConfusionPrecision) ) NbPoints=1;
2124
2125 //second side
2126 if (NbPoints<2) CalculPtsInterTriEdgeCoplanaires(TriSurfID,NormaleT,Tri1,Tri2,PE1,PE2,Edge,EdgeIndex,
2127 PT2,PT3,Cote23,2,SP1,SP2,NbPoints);
2128 }
2129 if ( (NbPoints>1)&&(Abs(SP2.U1()-SP1.U1())<MyConfusionPrecision)
2130 &&(Abs(SP2.V1()-SP1.V1())<MyConfusionPrecision) ) NbPoints=1;
2131 if (NbPoints>=2) return(NbPoints);
2132 }
2133 }
2134
2135 //------------------------------------------------------
2136 // NON COPLANAR edge and triangle (a contact point)
2137 //------------------------------------------------------
2138 else if( ( (pe1>=pt1)&&(pt1>=pe2) ) || ( (pe1<=pt1)&&(pt1<=pe2) ) ) { //
2139 lambda=(pe1-pt1)/(pe1-pe2);
2140 IntPolyh_Point PI;
2141 if (lambda<-MyConfusionPrecision) {
2142
2143 }
2144 else if (Abs(lambda)<MyConfusionPrecision) {//lambda==0
2145 PI=PE1;
2146 if(TriSurfID==1) SP1.SetEdge2(-1);
2147 else SP1.SetEdge1(-1);
2148 }
2149 else if (Abs(lambda-1.0)<MyConfusionPrecision) {//lambda==1
2150 PI=PE2;
2151 if(TriSurfID==1) SP1.SetEdge2(-1);
2152 else SP1.SetEdge1(-1);
2153 }
2154 else {
2155 PI=PE1+Edge*lambda;
2156 if(TriSurfID==1) {
2157 if(Tri2.GetEdgeOrientation(EdgeIndex)>0)
2158 SP1.SetLambda2(lambda);
2159 else SP1.SetLambda2(1.0-lambda);
2160 }
2161 if(TriSurfID==2) {
2162 if(Tri1.GetEdgeOrientation(EdgeIndex)>0)
2163 SP1.SetLambda1(lambda);
2164 else SP1.SetLambda1(1.0-lambda);
2165 }
2166 }
2167
2168 Standard_Real Cote23X=Cote23.X();
2169 Standard_Real D1=0.0;
2170 Standard_Real D3,D4;
2171
2172 //Combination Eq1 Eq2
2173 if(Abs(Cote23X)>MyConfusionPrecision) {
2174 D1=Cote12.Y()-Cote12.X()*Cote23.Y()/Cote23X;
2175 }
2176 if(Abs(D1)>MyConfusionPrecision) {
2177 alpha = ( PI.Y()-PT1.Y()-(PI.X()-PT1.X())*Cote23.Y()/Cote23X )/D1;
2178
2179 ///It is checked if 1.0>=alpha>=0.0
2180 if ((alpha<-MyConfusionPrecision)||(alpha>(1.0+MyConfusionPrecision))) return(0);
2181 else beta = (PI.X()-PT1.X()-alpha*Cote12.X())/Cote23X;
2182 }
2183 //Combination Eq1 and Eq2 with Cote23.X()==0
2184 else if ( (Abs(Cote12.X())>MyConfusionPrecision)
2185 &&(Abs(Cote23X)<MyConfusionPrecision) ) { //There is Cote23.X()==0
2186 alpha = (PI.X()-PT1.X())/Cote12.X();
2187
2188 if ((alpha<-MyConfusionPrecision)||(alpha>(1.0+MyConfusionPrecision))) return(0);
2189
2190 else if (Abs(Cote23.Y())>MyConfusionPrecision) beta = (PI.Y()-PT1.Y()-alpha*Cote12.Y())/Cote23.Y();
2191 else if (Abs(Cote23.Z())>MyConfusionPrecision) beta = (PI.Z()-PT1.Z()-alpha*Cote12.Z())/Cote23.Z();
2192 else {
2193
2194 }
2195 }
2196 //Combination Eq1 and Eq3
2197 else if ( (Abs(Cote23.X())>MyConfusionPrecision)
2198 &&(Abs( D3= (Cote12.Z()-Cote12.X()*Cote23.Z()/Cote23.X()) ) > MyConfusionPrecision) ) {
2199
2200 alpha = (PI.Z()-PT1.Z()-(PI.X()-PT1.X())*Cote23.Z()/Cote23.X())/D3;
2201
2202 if ( (alpha<-MyConfusionPrecision)||(alpha>(1.0+MyConfusionPrecision)) ) return(0);
2203 else beta = (PI.X()-PT1.X()-alpha*Cote12.X())/Cote23.X();
2204 }
2205 //Combination Eq2 and Eq3
2206 else if ( (Abs(Cote23.Y())>MyConfusionPrecision)
2207 &&(Abs( D4= (Cote12.Z()-Cote12.Y()*Cote23.Z()/Cote23.Y()) ) > MyConfusionPrecision) ) {
2208
2209 alpha = (PI.Z()-PT1.Z()-(PI.Y()-PT1.Y())*Cote23.Z()/Cote23.Y())/D4;
2210
2211 if ( (alpha<-MyConfusionPrecision)||(alpha>(1.0+MyConfusionPrecision)) ) return(0);
2212 else beta = (PI.Y()-PT1.Y()-alpha*Cote12.Y())/Cote23.Y();
2213 }
2214 //Combination Eq2 and Eq3 with Cote23.Y()==0
2215 else if ( (Abs(Cote12.Y())>MyConfusionPrecision)
2216 && (Abs(Cote23.Y())<MyConfusionPrecision) ) {
2217 alpha = (PI.Y()-PT1.Y())/Cote12.Y();
2218
2219 if ( (alpha<-MyConfusionPrecision)||(alpha>(1.0+MyConfusionPrecision)) ) return(0);
2220
2221 else if (Abs(Cote23.Z())>MyConfusionPrecision) beta = (PI.Z()-PT1.Z()-alpha*Cote12.Z())/Cote23.Z();
2222
2223 else {
2224 printf("\nCote PT2PT3 nul1\n");
2225 PT2.Dump(2004);
2226 PT3.Dump(3004);
2227 }
2228 }
2229 //Combination Eq1 and Eq3 with Cote23.Z()==0
2230 else if ( (Abs(Cote12.Z())>MyConfusionPrecision)
2231 && (Abs(Cote23.Z())<MyConfusionPrecision) ) {
2232 alpha = (PI.Z()-PT1.Z())/Cote12.Z();
2233
2234 if ( (alpha<-MyConfusionPrecision)||(alpha>(1.0+MyConfusionPrecision)) ) return(0);
2235
2236 else if (Abs(Cote23.X())>MyConfusionPrecision) beta = (PI.X()-PT1.X()-alpha*Cote12.X())/Cote23.X();
2237
2238 else {
2239
2240 }
2241 }
2242
2243 else { //Particular case not processed ?
2244
2245 alpha=RealLast();
2246 beta=RealLast();
2247 }
2248
2249 if( (beta<-MyConfusionPrecision)||(beta>(alpha+MyConfusionPrecision)) ) return(0);
2250 else {
2251 SP1.SetXYZ(PI.X(),PI.Y(),PI.Z());
2252
2253 if (TriSurfID==1) {
2254 SP1.SetUV2(PI.U(),PI.V());
2255 SP1.SetUV1(PT1.U()+Cote12.U()*alpha+Cote23.U()*beta, PT1.V()+Cote12.V()*alpha+Cote23.V()*beta);
2256 NbPoints++;
2257 if (alpha<MyConfusionPrecision) {//alpha=0 --> beta==0
2258 SP1.SetXYZ(PT1.X(),PT1.Y(),PT1.Z());
2259 SP1.SetUV1(PT1.U(),PT1.V());
2260 SP1.SetEdge1(-1);
2261 }
2262 else if ( (beta<MyConfusionPrecision)&&(Abs(1-alpha)<MyConfusionPrecision) ) {//beta==0 alpha==1
2263 SP1.SetXYZ(PT2.X(),PT2.Y(),PT2.Z());
2264 SP1.SetUV1(PT2.U(),PT2.V());
2265 SP1.SetEdge1(-1);
2266 }
2267 else if ( (Abs(beta-1)<MyConfusionPrecision)&&(Abs(1-alpha)<MyConfusionPrecision) ) {//beta==1 alpha==1
2268 SP1.SetXYZ(PT3.X(),PT3.Y(),PT3.Z());
2269 SP1.SetUV1(PT3.U(),PT3.V());
2270 SP1.SetEdge1(-1);
2271 }
2272 else if (beta<MyConfusionPrecision) {//beta==0
2273 SP1.SetEdge1(Tri1.GetEdgeNumber(1));
2274 if(Tri1.GetEdgeOrientation(1)>0)
2275 SP1.SetLambda1(alpha);
2276 else SP1.SetLambda1(1.0-alpha);
2277 }
2278 else if (Abs(beta-alpha)<MyConfusionPrecision) {//beta==alpha
2279 SP1.SetEdge1(Tri1.GetEdgeNumber(3));
2280 if(Tri1.GetEdgeOrientation(3)>0)
2281 SP1.SetLambda1(1.0-alpha);
2282 else SP1.SetLambda1(alpha);
2283 }
2284 else if (Abs(alpha-1)<MyConfusionPrecision) {//alpha==1
2285 SP1.SetEdge1(Tri1.GetEdgeNumber(2));
2286 if(Tri1.GetEdgeOrientation(2)>0)
2287 SP1.SetLambda1(beta);
2288 else SP1.SetLambda1(1.0-beta);
2289 }
2290 }
2291 else if(TriSurfID==2) {
2292 SP1.SetUV1(PI.U(),PI.V());
2293 SP1.SetUV2(PT1.U()+Cote12.U()*alpha+Cote23.U()*beta, PT1.V()+Cote12.V()*alpha+Cote23.V()*beta);
2294 NbPoints++;
2295 if (alpha<MyConfusionPrecision) {//alpha=0 --> beta==0
2296 SP1.SetXYZ(PT1.X(),PT1.Y(),PT1.Z());
2297 SP1.SetUV2(PT1.U(),PT1.V());
2298 SP1.SetEdge2(-1);
2299 }
2300 else if ( (beta<MyConfusionPrecision)&&(Abs(1-alpha)<MyConfusionPrecision) ) {//beta==0 alpha==1
2301 SP1.SetXYZ(PT2.X(),PT2.Y(),PT2.Z());
2302 SP1.SetUV2(PT2.U(),PT2.V());
2303 SP1.SetEdge2(-1);
2304 }
2305 else if ( (Abs(beta-1)<MyConfusionPrecision)&&(Abs(1-alpha)<MyConfusionPrecision) ) {//beta==1 alpha==1
2306 SP1.SetXYZ(PT3.X(),PT3.Y(),PT3.Z());
2307 SP1.SetUV2(PT3.U(),PT3.V());
2308 SP1.SetEdge2(-1);
2309 }
2310 else if (beta<MyConfusionPrecision) { //beta==0
2311 SP1.SetEdge2(Tri2.GetEdgeNumber(1));
2312 if(Tri2.GetEdgeOrientation(1)>0)
2313 SP1.SetLambda2(alpha);
2314 else SP1.SetLambda2(1.0-alpha);
2315 }
2316 else if (Abs(beta-alpha)<MyConfusionPrecision) {//beta==alpha
2317 SP1.SetEdge2(Tri2.GetEdgeNumber(3));
2318 if(Tri2.GetEdgeOrientation(3)>0)
2319 SP1.SetLambda2(1.0-alpha);
2320 else SP1.SetLambda2(alpha);
2321 }
2322 else if (Abs(alpha-1)<MyConfusionPrecision) {//alpha==1
2323 SP1.SetEdge2(Tri2.GetEdgeNumber(2));
2324 if(Tri2.GetEdgeOrientation(2)>0)
2325 SP1.SetLambda2(alpha);
2326 else SP1.SetLambda2(1.0-alpha);
2327 }
2328 }
2329 else {
2330
2331 }
2332 }
2333 }
2334 else return 0;
2335 }
2336 return (NbPoints);
2337 }
2338 //=======================================================================
2339 //function : TriangleCompare
2340 //purpose : Analyze each couple of triangles from the two --
2341 // array of triangles, to see if they are in
2342 // contact, and compute the incidence. Then put
2343 // couples in contact in the array of couples
2344 //=======================================================================
TriangleCompare()2345 Standard_Integer IntPolyh_MaillageAffinage::TriangleCompare ()
2346 {
2347 // Find couples with interfering bounding boxes
2348 IntPolyh_IndexedDataMapOfIntegerListOfInteger aDMILI;
2349 GetInterferingTriangles(TTriangles1, TPoints1,
2350 TTriangles2, TPoints2,
2351 aDMILI);
2352 if (aDMILI.IsEmpty()) {
2353 return 0;
2354 }
2355 //
2356 Standard_Real CoupleAngle = -2.0;
2357 //
2358 // Intersection of the triangles
2359 Standard_Integer i, aNb = aDMILI.Extent();
2360 for (i = 1; i <= aNb; ++i) {
2361 const Standard_Integer i_S1 = aDMILI.FindKey(i);
2362 IntPolyh_Triangle &Triangle1 = TTriangles1[i_S1];
2363 const IntPolyh_Point& P1 = TPoints1[Triangle1.FirstPoint()];
2364 const IntPolyh_Point& P2 = TPoints1[Triangle1.SecondPoint()];
2365 const IntPolyh_Point& P3 = TPoints1[Triangle1.ThirdPoint()];
2366 //
2367 const TColStd_ListOfInteger& aLI2 = aDMILI(i);
2368 TColStd_ListOfInteger::Iterator aItLI(aLI2);
2369 for (; aItLI.More(); aItLI.Next()) {
2370 const Standard_Integer i_S2 = aItLI.Value();
2371 IntPolyh_Triangle &Triangle2 = TTriangles2[i_S2];
2372 const IntPolyh_Point& Q1 = TPoints2[Triangle2.FirstPoint()];
2373 const IntPolyh_Point& Q2 = TPoints2[Triangle2.SecondPoint()];
2374 const IntPolyh_Point& Q3 = TPoints2[Triangle2.ThirdPoint()];
2375 //
2376 if (TriContact(P1, P2, P3, Q1, Q2, Q3, CoupleAngle)) {
2377 IntPolyh_Couple aCouple(i_S1, i_S2, CoupleAngle);
2378 TTrianglesContacts.Append(aCouple);
2379 //
2380 Triangle1.SetIntersection(Standard_True);
2381 Triangle2.SetIntersection(Standard_True);
2382 }
2383 }
2384 }
2385 return TTrianglesContacts.Extent();
2386 }
2387
2388 //=======================================================================
2389 //function : CheckCoupleAndGetAngle
2390 //purpose :
2391 //=======================================================================
CheckCoupleAndGetAngle(const Standard_Integer T1,const Standard_Integer T2,Standard_Real & Angle,IntPolyh_ListOfCouples & TTrianglesContacts)2392 Standard_Boolean CheckCoupleAndGetAngle(const Standard_Integer T1,
2393 const Standard_Integer T2,
2394 Standard_Real& Angle,
2395 IntPolyh_ListOfCouples &TTrianglesContacts)
2396 {
2397 IntPolyh_ListIteratorOfListOfCouples aIt(TTrianglesContacts);
2398 for (; aIt.More(); aIt.Next()) {
2399 IntPolyh_Couple& TestCouple = aIt.ChangeValue();
2400 if (!TestCouple.IsAnalyzed()) {
2401 if (TestCouple.FirstValue() == T1 && TestCouple.SecondValue() == T2) {
2402 TestCouple.SetAnalyzed(Standard_True);
2403 Angle = TestCouple.Angle();
2404 return Standard_True;
2405 }
2406 }
2407 }
2408 return Standard_False;
2409 }
2410 //=======================================================================
2411 //function : CheckCoupleAndGetAngle2
2412 //purpose :
2413 //=======================================================================
CheckCoupleAndGetAngle2(const Standard_Integer T1,const Standard_Integer T2,const Standard_Integer T11,const Standard_Integer T22,IntPolyh_ListIteratorOfListOfCouples & theItCT11,IntPolyh_ListIteratorOfListOfCouples & theItCT22,Standard_Real & Angle,IntPolyh_ListOfCouples & TTrianglesContacts)2414 Standard_Boolean CheckCoupleAndGetAngle2(const Standard_Integer T1,
2415 const Standard_Integer T2,
2416 const Standard_Integer T11,
2417 const Standard_Integer T22,
2418 IntPolyh_ListIteratorOfListOfCouples& theItCT11,
2419 IntPolyh_ListIteratorOfListOfCouples& theItCT22,
2420 Standard_Real & Angle,
2421 IntPolyh_ListOfCouples &TTrianglesContacts)
2422 {
2423 ///couple T1 T2 is found in the list
2424 ///T11 and T22 are two other triangles implied in the contact edge edge
2425 /// CT11 couple( T1,T22) and CT22 couple (T2,T11)
2426 /// these couples will be marked if there is a start point
2427 Standard_Boolean Test1 , Test2, Test3;
2428 Test1 = Test2 = Test3 = Standard_False;
2429 //
2430 IntPolyh_ListIteratorOfListOfCouples aIt(TTrianglesContacts);
2431 for (; aIt.More(); aIt.Next()) {
2432 IntPolyh_Couple& TestCouple = aIt.ChangeValue();
2433 if (TestCouple.IsAnalyzed()) {
2434 continue;
2435 }
2436 //
2437 if (TestCouple.FirstValue() == T1) {
2438 if (TestCouple.SecondValue() == T2) {
2439 Test1 = Standard_True;
2440 TestCouple.SetAnalyzed(Standard_True);
2441 Angle = TestCouple.Angle();
2442 }
2443 else if (TestCouple.SecondValue() == T22) {
2444 Test2 = Standard_True;
2445 theItCT11 = aIt;
2446 Angle = TestCouple.Angle();
2447 }
2448 }
2449 else if (TestCouple.FirstValue() == T11) {
2450 if (TestCouple.SecondValue() == T2) {
2451 Test3 = Standard_True;
2452 theItCT22 = aIt;
2453 Angle = TestCouple.Angle();
2454 }
2455 }
2456 //
2457 if (Test1 && Test2 && Test3) {
2458 break;
2459 }
2460 }
2461 return Test1;
2462 }
2463 //=======================================================================
2464 //function : CheckNextStartPoint
2465 //purpose : it is checked if the point is not a top
2466 // then it is stored in one or several valid arrays with
2467 // the proper list number
2468 //=======================================================================
CheckNextStartPoint(IntPolyh_SectionLine & SectionLine,IntPolyh_ArrayOfTangentZones & TTangentZones,IntPolyh_StartPoint & SP,const Standard_Boolean Prepend)2469 Standard_Integer CheckNextStartPoint(IntPolyh_SectionLine & SectionLine,
2470 IntPolyh_ArrayOfTangentZones & TTangentZones,
2471 IntPolyh_StartPoint & SP,
2472 const Standard_Boolean Prepend)//=Standard_False)
2473 {
2474 Standard_Integer Test=1;
2475 if( (SP.E1()==-1)||(SP.E2()==-1) ) {
2476 //The tops of triangle are analyzed
2477 //It is checked if they are not in the array TTangentZones
2478 Standard_Integer FinTTZ=TTangentZones.NbItems();
2479 for(Standard_Integer uiui=0; uiui<FinTTZ; uiui++) {
2480 IntPolyh_StartPoint TestSP=TTangentZones[uiui];
2481 if ( (Abs(SP.U1()-TestSP.U1())<MyConfusionPrecision)
2482 &&(Abs(SP.V1()-TestSP.V1())<MyConfusionPrecision) ) {
2483 if ( (Abs(SP.U2()-TestSP.U2())<MyConfusionPrecision)
2484 &&(Abs(SP.V2()-TestSP.V2())<MyConfusionPrecision) ) {
2485 Test=0;//SP is already in the list of tops
2486 uiui=FinTTZ;
2487 }
2488 }
2489 }
2490 if (Test) {//the top does not belong to the list of TangentZones
2491 SP.SetChainList(-1);
2492 TTangentZones[FinTTZ]=SP;
2493 TTangentZones.IncrementNbItems();
2494 Test=0;//the examined point is a top
2495 }
2496 }
2497 else if (Test) {
2498 if (Prepend)
2499 SectionLine.Prepend(SP);
2500 else {
2501 SectionLine[SectionLine.NbStartPoints()]=SP;
2502 SectionLine.IncrementNbStartPoints();
2503 }
2504
2505 }
2506 //if the point is not a top Test=1
2507 //The chain is continued
2508 return(Test);
2509 }
2510 //=======================================================================
2511 //function : StartPointsChain
2512 //purpose : Loop on the array of couples. Compute StartPoints.
2513 // Try to chain the StartPoints into SectionLines or
2514 // put the point in the ArrayOfTangentZones if
2515 // chaining it, is not possible.
2516 //=======================================================================
StartPointsChain(IntPolyh_ArrayOfSectionLines & TSectionLines,IntPolyh_ArrayOfTangentZones & TTangentZones)2517 Standard_Integer IntPolyh_MaillageAffinage::StartPointsChain
2518 (IntPolyh_ArrayOfSectionLines& TSectionLines,
2519 IntPolyh_ArrayOfTangentZones& TTangentZones)
2520 {
2521 //Loop on the array of couples filled in the function COMPARE()
2522 IntPolyh_ListIteratorOfListOfCouples aIt(TTrianglesContacts);
2523 for (; aIt.More(); aIt.Next()) {
2524 IntPolyh_Couple& aCouple = aIt.ChangeValue();
2525 // Check if the couple of triangles has not been already examined.
2526 if(!aCouple.IsAnalyzed()) {
2527
2528 Standard_Integer SectionLineIndex=TSectionLines.NbItems();
2529 // fill last section line if still empty (eap)
2530 if (SectionLineIndex > 0
2531 &&
2532 TSectionLines[SectionLineIndex-1].NbStartPoints() == 0)
2533 SectionLineIndex -= 1;
2534 else
2535 TSectionLines.IncrementNbItems();
2536
2537 IntPolyh_SectionLine & MySectionLine=TSectionLines[SectionLineIndex];
2538 if (MySectionLine.GetN() == 0) // eap
2539 MySectionLine.Init(10000);//Initialisation of array of StartPoint
2540
2541 Standard_Integer NbPoints=-1;
2542 Standard_Integer T1I, T2I;
2543 T1I = aCouple.FirstValue();
2544 T2I = aCouple.SecondValue();
2545
2546 // Start points for the current couple are found
2547 IntPolyh_StartPoint SP1, SP2;
2548 NbPoints=StartingPointsResearch(T1I,T2I,SP1, SP2);//first calculation
2549 aCouple.SetAnalyzed(Standard_True);//the couple is marked
2550
2551 if(NbPoints==1) {// particular case top/triangle or edge/edge
2552 //the start point is input in the array
2553 SP1.SetChainList(SectionLineIndex);
2554 SP1.SetAngle(aCouple.Angle());
2555 //it is checked if the point is not atop of the triangle
2556 if(CheckNextStartPoint(MySectionLine,TTangentZones,SP1)) {
2557 IntPolyh_StartPoint SPNext1;
2558 Standard_Integer TestSP1=0;
2559
2560 //chain of a side
2561 IntPolyh_StartPoint SP11;//=SP1;
2562 if(SP1.E1()>=0) { //&&(SP1.E2()!=-1) already tested if the point is not a top
2563 Standard_Integer NextTriangle1=0;
2564 if (TEdges1[SP1.E1()].FirstTriangle()!=T1I) NextTriangle1=TEdges1[SP1.E1()].FirstTriangle();
2565 else NextTriangle1=TEdges1[SP1.E1()].SecondTriangle();
2566
2567 Standard_Real Angle=-2.0;
2568 if (CheckCoupleAndGetAngle(NextTriangle1,T2I,Angle,TTrianglesContacts)) {
2569 //it is checked if the couple exists and is marked
2570 Standard_Integer NbPoints11=0;
2571 NbPoints11=NextStartingPointsResearch(NextTriangle1,T2I,SP1,SP11);
2572 if (NbPoints11==1) {
2573 SP11.SetChainList(SectionLineIndex);
2574 SP11.SetAngle(Angle);
2575
2576 if(CheckNextStartPoint(MySectionLine,TTangentZones,SP11)) {
2577 Standard_Integer EndChainList=1;
2578 while (EndChainList!=0) {
2579 TestSP1=GetNextChainStartPoint(SP11,SPNext1,MySectionLine,TTangentZones);
2580 if(TestSP1==1) {
2581 SPNext1.SetChainList(SectionLineIndex);
2582 if(CheckNextStartPoint(MySectionLine,TTangentZones,SPNext1))
2583 SP11=SPNext1;
2584 else EndChainList=0;
2585 }
2586 else EndChainList=0; //There is no next point
2587 }
2588 }
2589
2590 }
2591 else {
2592 if(NbPoints11>1) {//The point is input in the array TTangentZones
2593 TTangentZones[TTangentZones.NbItems()]=SP11;//default list number = -1
2594 TTangentZones.IncrementNbItems();
2595 }
2596 else {
2597
2598 }
2599 }
2600 }
2601 }
2602 else if (SP1.E2()<0){
2603
2604 }
2605 //chain of the other side
2606 IntPolyh_StartPoint SP12;//=SP1;
2607 if (SP1.E2()>=0) { //&&(SP1.E1()!=-1) already tested
2608 Standard_Integer NextTriangle2;
2609 if (TEdges2[SP1.E2()].FirstTriangle()!=T2I) NextTriangle2=TEdges2[SP1.E2()].FirstTriangle();
2610 else NextTriangle2=TEdges2[SP1.E2()].SecondTriangle();
2611
2612 Standard_Real Angle=-2.0;
2613 if(CheckCoupleAndGetAngle(T1I,NextTriangle2,Angle,TTrianglesContacts)) {
2614 Standard_Integer NbPoints12=0;
2615 NbPoints12=NextStartingPointsResearch(T1I,NextTriangle2,SP1, SP12);
2616 if (NbPoints12==1) {
2617
2618 SP12.SetChainList(SectionLineIndex);
2619 SP12.SetAngle(Angle);
2620 Standard_Boolean Prepend = Standard_True; // eap
2621
2622 if(CheckNextStartPoint(MySectionLine,TTangentZones,SP12, Prepend)) {
2623 Standard_Integer EndChainList=1;
2624 while (EndChainList!=0) {
2625 TestSP1=GetNextChainStartPoint(SP12,SPNext1,
2626 MySectionLine,TTangentZones,
2627 Prepend); // eap
2628 if(TestSP1==1) {
2629 SPNext1.SetChainList(SectionLineIndex);
2630 if(CheckNextStartPoint(MySectionLine,TTangentZones,SPNext1,Prepend))
2631 SP12=SPNext1;
2632 else EndChainList=0;
2633 }
2634 else EndChainList=0; //there is no next point
2635 }
2636 }
2637
2638 else {
2639 if(NbPoints12>1) {//The points are input in the array TTangentZones
2640 TTangentZones[TTangentZones.NbItems()]=SP12;//default list number = -1
2641 TTangentZones.IncrementNbItems();
2642 }
2643 else {
2644
2645 }
2646 }
2647 }
2648 }
2649 }
2650 else if(SP1.E1()<0){
2651
2652 }
2653 }
2654 }
2655 else if(NbPoints==2) {
2656 //the start points are input in the array
2657 IntPolyh_StartPoint SPNext2;
2658 Standard_Integer TestSP2=0;
2659 Standard_Integer EndChainList=1;
2660
2661 SP1.SetChainList(SectionLineIndex);
2662 SP1.SetAngle(aCouple.Angle());
2663 if(CheckNextStartPoint(MySectionLine,TTangentZones,SP1)) {
2664
2665 //chain of a side
2666 while (EndChainList!=0) {
2667 TestSP2=GetNextChainStartPoint(SP1,SPNext2,MySectionLine,TTangentZones);
2668 if(TestSP2==1) {
2669 SPNext2.SetChainList(SectionLineIndex);
2670 if(CheckNextStartPoint(MySectionLine,TTangentZones,SPNext2))
2671 SP1=SPNext2;
2672 else EndChainList=0;
2673 }
2674 else EndChainList=0; //there is no next point
2675 }
2676 }
2677
2678 SP2.SetChainList(SectionLineIndex);
2679 SP2.SetAngle(aCouple.Angle());
2680 Standard_Boolean Prepend = Standard_True; // eap
2681
2682 if(CheckNextStartPoint(MySectionLine,TTangentZones,SP2,Prepend)) {
2683
2684 //chain of the other side
2685 EndChainList=1;
2686 while (EndChainList!=0) {
2687 TestSP2=GetNextChainStartPoint(SP2,SPNext2,
2688 MySectionLine,TTangentZones,
2689 Prepend); // eap
2690 if(TestSP2==1) {
2691 SPNext2.SetChainList(SectionLineIndex);
2692 if(CheckNextStartPoint(MySectionLine,TTangentZones,SPNext2,Prepend))
2693 SP2=SPNext2;
2694 else EndChainList=0;
2695 }
2696 else EndChainList=0; //there is no next point
2697 }
2698 }
2699 }
2700
2701 else if( (NbPoints>2)&&(NbPoints<7) ) {
2702 //More than two start points
2703 //the start point is input in the table
2704 SP1.SetChainList(SectionLineIndex);
2705 CheckNextStartPoint(MySectionLine,TTangentZones,SP1);
2706 }
2707
2708 else {
2709
2710 }
2711 }
2712 }
2713
2714 return(1);
2715 }
2716 //=======================================================================
2717 //function : GetNextChainStartPoint
2718 //purpose : Mainly used by StartPointsChain(), this function
2719 // try to compute the next StartPoint.
2720 // GetNextChainStartPoint is used only if it is known that there are 2 contact points
2721 //=======================================================================
GetNextChainStartPoint(const IntPolyh_StartPoint & SP,IntPolyh_StartPoint & SPNext,IntPolyh_SectionLine & MySectionLine,IntPolyh_ArrayOfTangentZones & TTangentZones,const Standard_Boolean Prepend)2722 Standard_Integer IntPolyh_MaillageAffinage::GetNextChainStartPoint
2723 (const IntPolyh_StartPoint & SP,
2724 IntPolyh_StartPoint & SPNext,
2725 IntPolyh_SectionLine & MySectionLine,
2726 IntPolyh_ArrayOfTangentZones & TTangentZones,
2727 const Standard_Boolean Prepend)
2728 {
2729 Standard_Integer NbPoints=0;
2730 if( (SP.E1()>=0)&&(SP.E2()==-2) ) {
2731 //case if the point is on edge of T1
2732 Standard_Integer NextTriangle1;
2733 if (TEdges1[SP.E1()].FirstTriangle()!=SP.T1()) NextTriangle1=TEdges1[SP.E1()].FirstTriangle();
2734 else
2735 NextTriangle1=TEdges1[SP.E1()].SecondTriangle();
2736 //If is checked if two triangles intersect
2737 Standard_Real Angle= -2.0;
2738 if (CheckCoupleAndGetAngle(NextTriangle1,SP.T2(),Angle,TTrianglesContacts)) {
2739 NbPoints=NextStartingPointsResearch(NextTriangle1,SP.T2(),SP,SPNext);
2740 if( NbPoints!=1 ) {
2741 if (NbPoints>1)
2742 CheckNextStartPoint(MySectionLine,TTangentZones,SPNext,Prepend);
2743 else {
2744
2745 NbPoints=0;
2746 }
2747 }
2748 else
2749 SPNext.SetAngle(Angle);
2750 }
2751 else NbPoints=0;//this couple does not intersect
2752 }
2753 else if( (SP.E1()==-2)&&(SP.E2()>=0) ) {
2754 //case if the point is on edge of T2
2755 Standard_Integer NextTriangle2;
2756 if (TEdges2[SP.E2()].FirstTriangle()!=SP.T2()) NextTriangle2=TEdges2[SP.E2()].FirstTriangle();
2757 else
2758 NextTriangle2=TEdges2[SP.E2()].SecondTriangle();
2759 Standard_Real Angle= -2.0;
2760 if (CheckCoupleAndGetAngle(SP.T1(),NextTriangle2,Angle,TTrianglesContacts)) {
2761 NbPoints=NextStartingPointsResearch(SP.T1(),NextTriangle2,SP,SPNext);
2762 if( NbPoints!=1 ) {
2763 if (NbPoints>1)
2764 CheckNextStartPoint(MySectionLine,TTangentZones,SPNext,Prepend);
2765 else {
2766
2767 NbPoints=0;
2768 }
2769 }
2770 else
2771 SPNext.SetAngle(Angle);
2772 }
2773 else NbPoints=0;
2774 }
2775 else if( (SP.E1()==-2)&&(SP.E2()==-2) ) {
2776 ///no edge is touched or cut
2777
2778 NbPoints=0;
2779 }
2780 else if( (SP.E1()>=0)&&(SP.E2()>=0) ) {
2781 ///the point is located on two edges
2782 Standard_Integer NextTriangle1;
2783 if (TEdges1[SP.E1()].FirstTriangle()!=SP.T1()) NextTriangle1=TEdges1[SP.E1()].FirstTriangle();
2784 else
2785 NextTriangle1=TEdges1[SP.E1()].SecondTriangle();
2786 Standard_Integer NextTriangle2;
2787 if (TEdges2[SP.E2()].FirstTriangle()!=SP.T2()) NextTriangle2=TEdges2[SP.E2()].FirstTriangle();
2788 else
2789 NextTriangle2=TEdges2[SP.E2()].SecondTriangle();
2790 Standard_Real Angle= -2.0;
2791
2792 IntPolyh_ListIteratorOfListOfCouples aItCT11, aItCT22;
2793 if (CheckCoupleAndGetAngle2(NextTriangle1,NextTriangle2,
2794 SP.T1(),SP.T2(), aItCT11, aItCT22,
2795 Angle,TTrianglesContacts)) {
2796 NbPoints=NextStartingPointsResearch(NextTriangle1,NextTriangle2,SP,SPNext);
2797 if( NbPoints!=1 ) {
2798 if (NbPoints>1) {
2799 ///The new point is checked
2800 if(CheckNextStartPoint(MySectionLine,TTangentZones,SPNext,Prepend)>0) {
2801 }
2802 else {
2803
2804 }
2805 }
2806 NbPoints=0;
2807 }
2808 else {//NbPoints==1
2809 SPNext.SetAngle(Angle);
2810 if (aItCT11.More()) aItCT11.ChangeValue().SetAnalyzed(Standard_True);
2811 if (aItCT22.More()) aItCT22.ChangeValue().SetAnalyzed(Standard_True);
2812 }
2813 }
2814 else NbPoints=0;
2815 }
2816 else if( (SP.E1()==-1)||(SP.E2()==-1) ) {
2817 ///the points are tops of triangle
2818 ///the point is atored in an intermediary array
2819 }
2820 return(NbPoints);
2821 }
2822 //=======================================================================
2823 //function : GetArrayOfPoints
2824 //purpose :
2825 //=======================================================================
GetArrayOfPoints(const Standard_Integer SurfID) const2826 const IntPolyh_ArrayOfPoints& IntPolyh_MaillageAffinage::GetArrayOfPoints
2827 (const Standard_Integer SurfID)const
2828 {
2829 if (SurfID==1)
2830 return(TPoints1);
2831 return(TPoints2);
2832 }
2833 //=======================================================================
2834 //function : GetArrayOfEdges
2835 //purpose :
2836 //=======================================================================
GetArrayOfEdges(const Standard_Integer SurfID) const2837 const IntPolyh_ArrayOfEdges& IntPolyh_MaillageAffinage::GetArrayOfEdges
2838 (const Standard_Integer SurfID)const
2839 {
2840 if (SurfID==1)
2841 return(TEdges1);
2842 return(TEdges2);
2843 }
2844 //=======================================================================
2845 //function : GetArrayOfTriangles
2846 //purpose :
2847 //=======================================================================
2848 const IntPolyh_ArrayOfTriangles&
GetArrayOfTriangles(const Standard_Integer SurfID) const2849 IntPolyh_MaillageAffinage::GetArrayOfTriangles
2850 (const Standard_Integer SurfID)const{
2851 if (SurfID==1)
2852 return(TTriangles1);
2853 return(TTriangles2);
2854 }
2855
2856 //=======================================================================
2857 //function : GetBox
2858 //purpose :
2859 //=======================================================================
GetBox(const Standard_Integer SurfID) const2860 Bnd_Box IntPolyh_MaillageAffinage::GetBox(const Standard_Integer SurfID) const
2861 {
2862 if (SurfID==1)
2863 return(MyBox1);
2864 return(MyBox2);
2865 }
2866 //=======================================================================
2867 //function : GetArrayOfCouples
2868 //purpose :
2869 //=======================================================================
GetCouples()2870 IntPolyh_ListOfCouples &IntPolyh_MaillageAffinage::GetCouples()
2871 {
2872 return TTrianglesContacts;
2873 }
2874 //=======================================================================
2875 //function : SetEnlargeZone
2876 //purpose :
2877 //=======================================================================
SetEnlargeZone(const Standard_Boolean EnlargeZone)2878 void IntPolyh_MaillageAffinage::SetEnlargeZone(const Standard_Boolean EnlargeZone)
2879 {
2880 myEnlargeZone = EnlargeZone;
2881 }
2882 //=======================================================================
2883 //function : GetEnlargeZone
2884 //purpose :
2885 //=======================================================================
GetEnlargeZone() const2886 Standard_Boolean IntPolyh_MaillageAffinage::GetEnlargeZone() const
2887 {
2888 return myEnlargeZone;
2889 }
2890
2891 //=======================================================================
2892 //function : GetMinDeflection
2893 //purpose :
2894 //=======================================================================
GetMinDeflection(const Standard_Integer SurfID) const2895 Standard_Real IntPolyh_MaillageAffinage::GetMinDeflection(const Standard_Integer SurfID) const
2896 {
2897 return (SurfID==1)? FlecheMin1:FlecheMin2;
2898 }
2899
2900 //=======================================================================
2901 //function : GetMaxDeflection
2902 //purpose :
2903 //=======================================================================
GetMaxDeflection(const Standard_Integer SurfID) const2904 Standard_Real IntPolyh_MaillageAffinage::GetMaxDeflection(const Standard_Integer SurfID) const
2905 {
2906 return (SurfID==1)? FlecheMax1:FlecheMax2;
2907 }
2908
2909 //=======================================================================
2910 //function : DegeneratedIndex
2911 //purpose :
2912 //=======================================================================
DegeneratedIndex(const TColStd_Array1OfReal & aXpars,const Standard_Integer aNbX,const Handle (Adaptor3d_Surface)& aS,const Standard_Integer aIsoDirection,Standard_Integer & aI1,Standard_Integer & aI2)2913 void DegeneratedIndex(const TColStd_Array1OfReal& aXpars,
2914 const Standard_Integer aNbX,
2915 const Handle(Adaptor3d_Surface)& aS,
2916 const Standard_Integer aIsoDirection,
2917 Standard_Integer& aI1,
2918 Standard_Integer& aI2)
2919 {
2920 Standard_Integer i;
2921 Standard_Boolean bDegX1, bDegX2;
2922 Standard_Real aDegX1, aDegX2, aTol2, aX;
2923 //
2924 aI1=0;
2925 aI2=0;
2926 aTol2=MyTolerance*MyTolerance;
2927 //
2928 if (aIsoDirection==1){ // V=const
2929 bDegX1=IsDegenerated(aS, 1, aTol2, aDegX1);
2930 bDegX2=IsDegenerated(aS, 2, aTol2, aDegX2);
2931 }
2932 else if (aIsoDirection==2){ // U=const
2933 bDegX1=IsDegenerated(aS, 3, aTol2, aDegX1);
2934 bDegX2=IsDegenerated(aS, 4, aTol2, aDegX2);
2935 }
2936 else {
2937 return;
2938 }
2939 //
2940 if (!(bDegX1 || bDegX2)) {
2941 return;
2942 }
2943 //
2944 for(i=1; i<=aNbX; ++i) {
2945 aX=aXpars(i);
2946 if (bDegX1) {
2947 if (fabs(aX-aDegX1) < MyTolerance) {
2948 aI1=i;
2949 }
2950 }
2951 if (bDegX2) {
2952 if (fabs(aX-aDegX2) < MyTolerance) {
2953 aI2=i;
2954 }
2955 }
2956 }
2957 }
2958 //=======================================================================
2959 //function : IsDegenerated
2960 //purpose :
2961 //=======================================================================
IsDegenerated(const Handle (Adaptor3d_Surface)& aS,const Standard_Integer aIndex,const Standard_Real aTol2,Standard_Real & aDegX)2962 Standard_Boolean IsDegenerated(const Handle(Adaptor3d_Surface)& aS,
2963 const Standard_Integer aIndex,
2964 const Standard_Real aTol2,
2965 Standard_Real& aDegX)
2966 {
2967 Standard_Boolean bRet;
2968 Standard_Integer i, aNbP;
2969 Standard_Real aU, dU, aU1, aU2, aV, dV, aV1, aV2, aD2;
2970 gp_Pnt aP1, aP2;
2971 //
2972 bRet=Standard_False;
2973 aNbP=3;
2974 aDegX=99;
2975 //
2976 aU1=aS->FirstUParameter();
2977 aU2=aS->LastUParameter();
2978 aV1=aS->FirstVParameter();
2979 aV2=aS->LastVParameter();
2980 //
2981 if (aIndex<3) { // V=const
2982 aV=aV1;
2983 if (aIndex==2) {
2984 aV=aV2;
2985 }
2986 dU=(aU2-aU1)/(aNbP-1);
2987 aU=aU1;
2988 aP1=aS->Value(aU, aV);
2989 for (i=1; i<aNbP; ++i) {
2990 aU=i*dU;
2991 if (i==aNbP-1){
2992 aU=aU2;
2993 }
2994 aP2=aS->Value(aU, aV);
2995 aD2=aP1.SquareDistance(aP2);
2996 if (aD2>aTol2) {
2997 return bRet;
2998 }
2999 aP1=aP2;
3000 }
3001 aDegX=aV;
3002 bRet=!bRet;
3003 }
3004 else {// U=const
3005 aU=aU1;
3006 if (aIndex==4) {
3007 aU=aU2;
3008 }
3009 dV=(aV2-aV1)/(aNbP-1);
3010 aV=aV1;
3011 aP1=aS->Value(aU, aV);
3012 for (i=1; i<aNbP; ++i) {
3013 aV=i*dV;
3014 if (i==aNbP-1){
3015 aV=aV2;
3016 }
3017 aP2=aS->Value(aU, aV);
3018 aD2=aP1.SquareDistance(aP2);
3019 if (aD2>aTol2) {
3020 return bRet;
3021 }
3022 aP1=aP2;
3023 }
3024 bRet=!bRet;
3025 aDegX=aU;
3026 }
3027 //
3028 return bRet;
3029 }
3030