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
6 // This file is part of Open CASCADE Technology software library.
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.
14 // Alternatively, this file may be used under the terms of Open CASCADE
15 // commercial license or contractual agreement.
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
26 #include <Adaptor3d_HSurface.hxx>
27 #include <Bnd_BoundSortBox.hxx>
28 #include <Bnd_Box.hxx>
29 #include <Bnd_HArray1OfBox.hxx>
32 #include <IntCurveSurface_ThePolyhedronOfHInter.hxx>
33 #include <IntPolyh_ArrayOfCouples.hxx>
34 #include <IntPolyh_Couple.hxx>
35 #include <IntPolyh_Edge.hxx>
36 #include <IntPolyh_MaillageAffinage.hxx>
37 #include <IntPolyh_Point.hxx>
38 #include <IntPolyh_SectionLine.hxx>
39 #include <IntPolyh_StartPoint.hxx>
40 #include <IntPolyh_Triangle.hxx>
41 #include <Precision.hxx>
42 #include <TColStd_ListIteratorOfListOfInteger.hxx>
44 static Standard_Real MyTolerance=10.0e-7;
45 static Standard_Real MyConfusionPrecision=10.0e-12;
46 static Standard_Real SquareMyConfusionPrecision=10.0e-24;
49 inline Standard_Real maxSR(const Standard_Real a,
50 const Standard_Real b,
51 const Standard_Real c);
54 inline Standard_Real minSR(const Standard_Real a,
55 const Standard_Real b,
56 const Standard_Real c);
58 Standard_Integer project6(const IntPolyh_Point &ax,
59 const IntPolyh_Point &p1,
60 const IntPolyh_Point &p2,
61 const IntPolyh_Point &p3,
62 const IntPolyh_Point &q1,
63 const IntPolyh_Point &q2,
64 const IntPolyh_Point &q3);
66 void TestNbPoints(const Standard_Integer ,
67 Standard_Integer &NbPoints,
68 Standard_Integer &NbPointsTotal,
69 const IntPolyh_StartPoint &Pt1,
70 const IntPolyh_StartPoint &Pt2,
71 IntPolyh_StartPoint &SP1,
72 IntPolyh_StartPoint &SP2);
74 void CalculPtsInterTriEdgeCoplanaires(const Standard_Integer TriSurfID,
75 const IntPolyh_Point &NormaleTri,
76 const IntPolyh_Point &PE1,
77 const IntPolyh_Point &PE2,
78 const IntPolyh_Point &Edge,
79 const IntPolyh_Point &PT1,
80 const IntPolyh_Point &PT2,
81 const IntPolyh_Point &Cote,
82 const Standard_Integer CoteIndex,
83 IntPolyh_StartPoint &SP1,
84 IntPolyh_StartPoint &SP2,
85 Standard_Integer &NbPoints);
87 void CalculPtsInterTriEdgeCoplanaires2(const Standard_Integer TriSurfID,
88 const IntPolyh_Point &NormaleTri,
89 const IntPolyh_Triangle &Tri1,
90 const IntPolyh_Triangle &Tri2,
91 const IntPolyh_Point &PE1,
92 const IntPolyh_Point &PE2,
93 const IntPolyh_Point &Edge,
94 const Standard_Integer EdgeIndex,
95 const IntPolyh_Point &PT1,
96 const IntPolyh_Point &PT2,
97 const IntPolyh_Point &Cote,
98 const Standard_Integer CoteIndex,
99 IntPolyh_StartPoint &SP1,
100 IntPolyh_StartPoint &SP2,
101 Standard_Integer &NbPoints);
103 Standard_Integer CheckCoupleAndGetAngle(const Standard_Integer T1,
104 const Standard_Integer T2,
105 Standard_Real& Angle,
106 IntPolyh_ArrayOfCouples &TTrianglesContacts);
108 Standard_Integer CheckCoupleAndGetAngle2(const Standard_Integer T1,
109 const Standard_Integer T2,
110 const Standard_Integer T11,
111 const Standard_Integer T22,
112 Standard_Integer &CT11,
113 Standard_Integer &CT22,
114 Standard_Real & Angle,
115 IntPolyh_ArrayOfCouples &TTrianglesContacts);
117 Standard_Integer CheckNextStartPoint(IntPolyh_SectionLine & SectionLine,
118 IntPolyh_ArrayOfTangentZones & TTangentZones,
119 IntPolyh_StartPoint & SP,
120 const Standard_Boolean Prepend=Standard_False);
123 Standard_Boolean IsDegenerated(const Handle(Adaptor3d_HSurface)& aS,
124 const Standard_Integer aIndex,
125 const Standard_Real aTol2,
126 Standard_Real& aDegX);
128 void DegeneratedIndex(const TColStd_Array1OfReal& Xpars,
129 const Standard_Integer aNbX,
130 const Handle(Adaptor3d_HSurface)& aS,
131 const Standard_Integer aIsoDirection,
132 Standard_Integer& aI1,
133 Standard_Integer& aI2);
135 void EnlargeZone(const Handle(Adaptor3d_HSurface)& MaSurface,
141 //=======================================================================
142 //function : IntPolyh_MaillageAffinage
144 //=======================================================================
145 IntPolyh_MaillageAffinage::IntPolyh_MaillageAffinage
146 (const Handle(Adaptor3d_HSurface)& Surface1,
147 const Handle(Adaptor3d_HSurface)& Surface2,
148 const Standard_Integer )
150 MaSurface1(Surface1),
151 MaSurface2(Surface2),
162 myEnlargeZone(Standard_False)
165 //=======================================================================
166 //function : IntPolyh_MaillageAffinage
168 //=======================================================================
169 IntPolyh_MaillageAffinage::IntPolyh_MaillageAffinage
170 (const Handle(Adaptor3d_HSurface)& Surface1,
171 const Standard_Integer NbSU1,
172 const Standard_Integer NbSV1,
173 const Handle(Adaptor3d_HSurface)& Surface2,
174 const Standard_Integer NbSU2,
175 const Standard_Integer NbSV2,
176 const Standard_Integer )
178 MaSurface1(Surface1),
179 MaSurface2(Surface2),
190 myEnlargeZone(Standard_False)
193 //=======================================================================
194 //function : FillArrayOfPnt
195 //purpose : Compute points on one surface and fill an array of points
196 //=======================================================================
197 void IntPolyh_MaillageAffinage::FillArrayOfPnt
198 (const Standard_Integer SurfID)
200 Standard_Integer NbSamplesU, NbSamplesV, i, aNbSamplesU1, aNbSamplesV1;
201 Standard_Real u0, u1, v0, v1, aU, aV, dU, dV;
203 const Handle(Adaptor3d_HSurface)& MaSurface=(SurfID==1)? MaSurface1 : MaSurface2;
204 NbSamplesU=(SurfID==1)? NbSamplesU1:NbSamplesU2;
205 NbSamplesV=(SurfID==1)? NbSamplesV1:NbSamplesV2;
207 u0 = (MaSurface)->FirstUParameter();
208 u1 = (MaSurface)->LastUParameter();
209 v0 = (MaSurface)->FirstVParameter();
210 v1 = (MaSurface)->LastVParameter();
213 EnlargeZone(MaSurface, u0, u1, v0, v1);
216 TColStd_Array1OfReal aUpars(1, NbSamplesU);
217 TColStd_Array1OfReal aVpars(1, NbSamplesV);
219 aNbSamplesU1=NbSamplesU-1;
220 aNbSamplesV1=NbSamplesV-1;
222 dU=(u1-u0)/Standard_Real(aNbSamplesU1);
223 dV=(v1-v0)/Standard_Real(aNbSamplesV1);
225 for (i=0; i<NbSamplesU; ++i) {
227 if (i==aNbSamplesU1) {
230 aUpars.SetValue(i+1, aU);
233 for (i=0; i<NbSamplesV; ++i) {
235 if (i==aNbSamplesV1) {
238 aVpars.SetValue(i+1, aV);
241 FillArrayOfPnt(SurfID, aUpars, aVpars);
243 //=======================================================================
244 //function : FillArrayOfPnt
245 //purpose : Compute points on one surface and fill an array of points
246 // FILL AN ARRAY OF POINTS
247 //=======================================================================
248 void IntPolyh_MaillageAffinage::FillArrayOfPnt
249 (const Standard_Integer SurfID,
250 const Standard_Boolean isShiftFwd)
252 Standard_Integer NbSamplesU, NbSamplesV, i, aNbSamplesU1, aNbSamplesV1;
253 Standard_Real u0, u1, v0, v1, aU, aV, dU, dV;
254 const Handle(Adaptor3d_HSurface)& MaSurface=(SurfID==1)? MaSurface1 : MaSurface2;
255 NbSamplesU=(SurfID==1)? NbSamplesU1:NbSamplesU2;
256 NbSamplesV=(SurfID==1)? NbSamplesV1:NbSamplesV2;
258 u0 = (MaSurface)->FirstUParameter();
259 u1 = (MaSurface)->LastUParameter();
260 v0 = (MaSurface)->FirstVParameter();
261 v1 = (MaSurface)->LastVParameter();
264 EnlargeZone(MaSurface, u0, u1, v0, v1);
267 TColStd_Array1OfReal aUpars(1, NbSamplesU);
268 TColStd_Array1OfReal aVpars(1, NbSamplesV);
270 aNbSamplesU1=NbSamplesU-1;
271 aNbSamplesV1=NbSamplesV-1;
273 dU=(u1-u0)/Standard_Real(aNbSamplesU1);
274 dV=(v1-v0)/Standard_Real(aNbSamplesV1);
276 for (i=0; i<NbSamplesU; ++i) {
278 if (i==aNbSamplesU1) {
281 aUpars.SetValue(i+1, aU);
284 for (i=0; i<NbSamplesV; ++i) {
286 if (i==aNbSamplesV1) {
289 aVpars.SetValue(i+1, aV);
292 FillArrayOfPnt(SurfID, isShiftFwd, aUpars, aVpars);
294 //=======================================================================
295 //function : FillArrayOfPnt
296 //purpose : Compute points on one surface and fill an array of points
297 //=======================================================================
298 void IntPolyh_MaillageAffinage::FillArrayOfPnt
299 (const Standard_Integer SurfID,
300 const TColStd_Array1OfReal& Upars,
301 const TColStd_Array1OfReal& Vpars)
303 Standard_Boolean bDegI, bDeg;
304 Standard_Integer aNbU, aNbV, iCnt, i, j;
305 Standard_Integer aID1, aID2, aJD1, aJD2;
306 Standard_Real aTol, aU, aV, aX, aY, aZ;
309 aNbU=(SurfID==1)? NbSamplesU1 : NbSamplesU2;
310 aNbV=(SurfID==1)? NbSamplesV1 : NbSamplesV2;
311 Bnd_Box& aBox = (SurfID==1) ? MyBox1 : MyBox2;
312 Handle(Adaptor3d_HSurface)& aS=(SurfID==1)? MaSurface1:MaSurface2;
313 IntPolyh_ArrayOfPoints &TPoints=(SurfID==1)? TPoints1:TPoints2;
319 DegeneratedIndex(Vpars, aNbV, aS, 1, aJD1, aJD2);
320 if (!(aJD1 || aJD2)) {
321 DegeneratedIndex(Upars, aNbU, aS, 2, aID1, aID2);
324 TPoints.Init(aNbU*aNbV);
326 for(i=1; i<=aNbU; ++i){
327 bDegI=(aID1==i || aID2==i);
329 for(j=1; j<=aNbV; ++j){
331 aP=aS->Value(aU, aV);
332 aP.Coord(aX, aY, aZ);
333 IntPolyh_Point& aIP=TPoints[iCnt];
334 aIP.Set(aX, aY, aZ, aU, aV);
336 bDeg=bDegI || (aJD1==j || aJD2==j);
338 aIP.SetDegenerated(bDeg);
345 TPoints.SetNbItems(iCnt);
347 IntCurveSurface_ThePolyhedronOfHInter polyhedron(aS, Upars, Vpars);
349 aTol=polyhedron.DeflectionOverEstimation();
352 Standard_Real a1,a2,a3,b1,b2,b3;
354 aBox.Get(a1,a2,a3,b1,b2,b3);
355 aBox.Update(a1-aTol,a2-aTol,a3-aTol,b1+aTol,b2+aTol,b3+aTol);
356 aBox.Enlarge(MyTolerance);
359 //=======================================================================
360 //function : FillArrayOfPnt
361 //purpose : Compute points on one surface and fill an array of points
362 // REMPLISSAGE DU TABLEAU DE POINTS
363 //=======================================================================
364 void IntPolyh_MaillageAffinage::FillArrayOfPnt
365 (const Standard_Integer SurfID,
366 const Standard_Boolean isShiftFwd,
367 const TColStd_Array1OfReal& Upars,
368 const TColStd_Array1OfReal& Vpars)
370 Standard_Boolean bDegI, bDeg;
371 Standard_Integer aNbU, aNbV, iCnt, i, j;
372 Standard_Integer aID1, aID2, aJD1, aJD2;
373 Standard_Real Tol, resol, aU, aV, aMag;
374 Standard_Real aX, aY, aZ;
376 gp_Vec aDU, aDV, aNorm;
378 aNbU=(SurfID==1)? NbSamplesU1:NbSamplesU2;
379 aNbV=(SurfID==1)? NbSamplesV1:NbSamplesV2;
380 Bnd_Box& aBox = (SurfID==1) ? MyBox1 : MyBox2;
381 Handle(Adaptor3d_HSurface) aS=(SurfID==1)? MaSurface1:MaSurface2;
382 IntPolyh_ArrayOfPoints &TPoints=(SurfID==1)? TPoints1:TPoints2;
384 resol = gp::Resolution();
386 IntCurveSurface_ThePolyhedronOfHInter polyhedron(aS, Upars, Vpars);
387 Tol=polyhedron.DeflectionOverEstimation();
392 DegeneratedIndex(Vpars, aNbV, aS, 1, aJD1, aJD2);
393 if (!(aJD1 || aJD2)) {
394 DegeneratedIndex(Upars, aNbU, aS, 2, aID1, aID2);
397 TPoints.Init(aNbU*aNbV);
399 for(i=1; i<=aNbU; ++i){
400 bDegI=(aID1==i || aID2==i);
402 for(j=1; j<=aNbV; ++j){
404 aS->D1(aU, aV, aP, aDU, aDV);
406 aNorm = aDU.Crossed(aDV);
407 aMag = aNorm.Magnitude();
410 aNorm.Multiply(Tol*1.5);
416 aP.Translate(aNorm.Reversed());
420 IntPolyh_Point& aIP=TPoints[iCnt];
421 aP.Coord(aX, aY, aZ);
422 aIP.Set(aX, aY, aZ, aU, aV);
424 bDeg=bDegI || (aJD1==j || aJD2==j);
426 aIP.SetDegenerated(bDeg);
433 TPoints.SetNbItems(iCnt);
437 Standard_Real a1,a2,a3,b1,b2,b3;
439 aBox.Get(a1,a2,a3,b1,b2,b3);
440 aBox.Update(a1-Tol,a2-Tol,a3-Tol,b1+Tol,b2+Tol,b3+Tol);
441 aBox.Enlarge(MyTolerance);
443 //=======================================================================
444 //function : CommonBox
445 //purpose : Compute the common box witch is the intersection
446 // of the two bounding boxes, and mark the points of
447 // the two surfaces that are inside.
448 // REJECTION BOUNDING BOXES
449 // DETERMINATION OF THE COMMON BOX
450 //=======================================================================
451 void IntPolyh_MaillageAffinage::CommonBox (const Bnd_Box &,
460 Standard_Real x10,y10,z10,x11,y11,z11;
461 Standard_Real x20,y20,z20,x21,y21,z21;
463 MyBox1.Get(x10,y10,z10,x11,y11,z11);
464 MyBox2.Get(x20,y20,z20,x21,y21,z21);
472 if((x10>x21)||(x20>x11)||(y10>y21)||(y20>y11)||(z10>z21)||(z20>z11)) {
475 if(x11<=x21) XMax=x11; else { if(x21<=x11) XMax=x21;}
476 if(x20<=x10) XMin=x10; else { if(x10<=x20) XMin=x20;}
477 if(y11<=y21) YMax=y11; else { if(y21<=y11) YMax=y21;}
478 if(y20<=y10) YMin=y10; else { if(y10<=y20) YMin=y20;}
479 if(z11<=z21) ZMax=z11; else { if(z21<=z11) ZMax=z21;}
480 if(z20<=z10) ZMin=z10; else { if(z10<=z20) ZMin=z20;}
482 if(((XMin==XMax)&&(!(YMin==YMax)&&!(ZMin==ZMax)))
483 ||((YMin==YMax)&&(!(XMin==XMax)&&!(ZMin==ZMax)))//ou exclusif ??
484 ||((ZMin==ZMax)&&(!(XMin==XMax)&&!(YMin==YMax)))) {
492 //extension of the box
493 if( (X==0)&&(Y!=0) ) X=Y*0.1;
494 else if( (X==0)&&(Z!=0) ) X=Z*0.1;
497 if( (Y==0)&&(X!=0) ) Y=X*0.1;
498 else if( (Y==0)&&(Z!=0) ) Y=Z*0.1;
501 if( (Z==0)&&(X!=0) ) Z=X*0.1;
502 else if( (Z==0)&&(Y!=0) ) Z=Y*0.1;
506 if( (X==0)&&(Y==0)&&(Z==0) ) {
514 //Marking of points included in the common
515 const Standard_Integer FinTP1 = TPoints1.NbItems();
516 // for(Standard_Integer i=0; i<FinTP1; i++) {
518 for( i=0; i<FinTP1; i++) {
519 IntPolyh_Point & Pt1 = TPoints1[i];
546 Pt1.SetPartOfCommon(r);
549 const Standard_Integer FinTP2 = TPoints2.NbItems();
550 for(Standard_Integer ii=0; ii<FinTP2; ii++) {
551 IntPolyh_Point & Pt2 = TPoints2[ii];
579 Pt2.SetPartOfCommon(rr);
582 //=======================================================================
583 //function : FillArrayOfEdges
584 //purpose : Compute edges from the array of points
585 // FILL THE ARRAY OF EDGES
586 //=======================================================================
587 void IntPolyh_MaillageAffinage::FillArrayOfEdges
588 (const Standard_Integer SurfID)
591 IntPolyh_ArrayOfEdges &TEdges=(SurfID==1)? TEdges1:TEdges2;
592 Standard_Integer NbSamplesU=(SurfID==1)? NbSamplesU1:NbSamplesU2;
593 Standard_Integer NbSamplesV=(SurfID==1)? NbSamplesV1:NbSamplesV2;
595 //NbEdges = 3 + 3*(NbSamplesV-2) + 3*(NbSamplesU-2) +
596 // + 3*(NbSamplesU-2)*(NbSamplesV-2) + (NbSamplesV-1) + (NbSamplesU-1);
597 //NbSamplesU and NbSamples cannot be less than 2, so
598 Standard_Integer NbEdges = 3*NbSamplesU*NbSamplesV - 2*(NbSamplesU+NbSamplesV) + 1;
599 TEdges.Init(NbEdges);
601 Standard_Integer CpteurTabEdges=0;
604 TEdges[CpteurTabEdges].SetFirstPoint(0); // U V
605 TEdges[CpteurTabEdges].SetSecondPoint(1); // U V+1
606 // TEdges[CpteurTabEdges].SetFirstTriangle(-1);
607 TEdges[CpteurTabEdges].SetSecondTriangle(0);
610 TEdges[CpteurTabEdges].SetFirstPoint(0); // U V
611 TEdges[CpteurTabEdges].SetSecondPoint(NbSamplesV); // U+1 V
612 TEdges[CpteurTabEdges].SetFirstTriangle(0);
613 TEdges[CpteurTabEdges].SetSecondTriangle(1);
616 TEdges[CpteurTabEdges].SetFirstPoint(0); // U V
617 TEdges[CpteurTabEdges].SetSecondPoint(NbSamplesV+1); // U+1 V+1
618 TEdges[CpteurTabEdges].SetFirstTriangle(1);
619 // TEdges[CpteurTabEdges].SetSecondTriangle(-1);
623 Standard_Integer PntInit=1;
624 Standard_Integer BoucleMeshV;
625 for(BoucleMeshV=1; BoucleMeshV<NbSamplesV-1;BoucleMeshV++){
626 TEdges[CpteurTabEdges].SetFirstPoint(PntInit); // U V
627 TEdges[CpteurTabEdges].SetSecondPoint(PntInit+1); // U V+1
628 // TEdges[CpteurTabEdges].SetFirstTriangle(-1);
629 TEdges[CpteurTabEdges].SetSecondTriangle(BoucleMeshV*2);
632 TEdges[CpteurTabEdges].SetFirstPoint(PntInit); // U V
633 TEdges[CpteurTabEdges].SetSecondPoint(PntInit+NbSamplesV+1); // U+1 V+1
634 TEdges[CpteurTabEdges].SetFirstTriangle(BoucleMeshV*2);
635 TEdges[CpteurTabEdges].SetSecondTriangle(BoucleMeshV*2+1);
638 TEdges[CpteurTabEdges].SetFirstPoint(PntInit); // U V
639 TEdges[CpteurTabEdges].SetSecondPoint(PntInit+NbSamplesV); // U+1 V
640 TEdges[CpteurTabEdges].SetFirstTriangle(BoucleMeshV*2+1);
641 TEdges[CpteurTabEdges].SetSecondTriangle(BoucleMeshV*2-2);
648 for(BoucleMeshV=1; BoucleMeshV<NbSamplesU-1;BoucleMeshV++){
649 TEdges[CpteurTabEdges].SetFirstPoint(PntInit); // U V
650 TEdges[CpteurTabEdges].SetSecondPoint(PntInit+1); // U V+1
651 TEdges[CpteurTabEdges].SetFirstTriangle((BoucleMeshV-1)*(NbSamplesV-1)*2+1);
652 TEdges[CpteurTabEdges].SetSecondTriangle(BoucleMeshV*(NbSamplesV-1)*2);
655 TEdges[CpteurTabEdges].SetFirstPoint(PntInit); // U V
656 TEdges[CpteurTabEdges].SetSecondPoint(PntInit+NbSamplesV+1); // U+1 V+1
657 TEdges[CpteurTabEdges].SetFirstTriangle(BoucleMeshV*(NbSamplesV-1)*2);
658 TEdges[CpteurTabEdges].SetSecondTriangle(BoucleMeshV*(NbSamplesV-1)*2+1);
661 TEdges[CpteurTabEdges].SetFirstPoint(PntInit); // U V
662 TEdges[CpteurTabEdges].SetSecondPoint(PntInit+NbSamplesV); // U+1 V
663 TEdges[CpteurTabEdges].SetFirstTriangle(BoucleMeshV*(NbSamplesV-1)*2+1);
664 // TEdges[CpteurTabEdges].SetSecondTriangle(-1);
669 PntInit=NbSamplesV+1;
670 //To provide recursion I associate a point with three edges
671 for(Standard_Integer BoucleMeshU=1; BoucleMeshU<NbSamplesU-1; BoucleMeshU++){
672 for(BoucleMeshV=1; BoucleMeshV<NbSamplesV-1;BoucleMeshV++){
673 TEdges[CpteurTabEdges].SetFirstPoint(PntInit); // U V
674 TEdges[CpteurTabEdges].SetSecondPoint(PntInit+1); // U V+1
675 TEdges[CpteurTabEdges].SetFirstTriangle((NbSamplesV-1)*2*(BoucleMeshU-1)+BoucleMeshV*2+1);
676 TEdges[CpteurTabEdges].SetSecondTriangle((NbSamplesV-1)*2*BoucleMeshU+BoucleMeshV*2);
679 TEdges[CpteurTabEdges].SetFirstPoint(PntInit); // U V
680 TEdges[CpteurTabEdges].SetSecondPoint(PntInit+NbSamplesV+1); // U+1 V+1
681 TEdges[CpteurTabEdges].SetFirstTriangle((NbSamplesV-1)*2*BoucleMeshU+BoucleMeshV*2);
682 TEdges[CpteurTabEdges].SetSecondTriangle((NbSamplesV-1)*2*BoucleMeshU+BoucleMeshV*2+1);
685 TEdges[CpteurTabEdges].SetFirstPoint(PntInit); // U V
686 TEdges[CpteurTabEdges].SetSecondPoint(PntInit+NbSamplesV); // U+1 V
687 TEdges[CpteurTabEdges].SetFirstTriangle((NbSamplesV-1)*2*BoucleMeshU+BoucleMeshV*2+1);
688 TEdges[CpteurTabEdges].SetSecondTriangle((NbSamplesV-1)*2*BoucleMeshU+BoucleMeshV*2-2);
690 PntInit++;//Pass to the next point
692 PntInit++;//Pass the last point of the column
693 PntInit++;//Pass the first point of the next column
697 PntInit=(NbSamplesU-1)*NbSamplesV; //point U=u1 V=0
698 for(BoucleMeshV=0; BoucleMeshV<NbSamplesV-1; BoucleMeshV++){
699 TEdges[CpteurTabEdges].SetFirstPoint(PntInit); //U=u1 V
700 TEdges[CpteurTabEdges].SetSecondPoint(PntInit+1); //U=u1 V+1
701 TEdges[CpteurTabEdges].SetFirstTriangle((NbSamplesU-2)*(NbSamplesV-1)*2+BoucleMeshV*2+1);
702 // TEdges[CpteurTabEdges].SetSecondTriangle(-1);
708 for(BoucleMeshV=0; BoucleMeshV<NbSamplesU-1;BoucleMeshV++){
709 TEdges[CpteurTabEdges].SetFirstPoint(NbSamplesV-1+BoucleMeshV*NbSamplesV); // U V=v1
710 TEdges[CpteurTabEdges].SetSecondPoint(NbSamplesV-1+(BoucleMeshV+1)*NbSamplesV); //U+1 V=v1
711 // TEdges[CpteurTabEdges].SetFirstTriangle(-1);
712 TEdges[CpteurTabEdges].SetSecondTriangle(BoucleMeshV*2*(NbSamplesV-1)+(NbSamplesV-2)*2);
715 TEdges.SetNbItems(CpteurTabEdges);
718 //=======================================================================
719 //function : FillArrayOfTriangles
720 //purpose : Compute triangles from the array of points, and --
721 // mark the triangles that use marked points by the
722 // CommonBox function.
723 // FILL THE ARRAY OF TRIANGLES
724 //=======================================================================
725 void IntPolyh_MaillageAffinage::FillArrayOfTriangles
726 (const Standard_Integer SurfID)
728 Standard_Integer CpteurTabTriangles=0;
729 Standard_Integer PntInit=0;
731 IntPolyh_ArrayOfPoints &TPoints=(SurfID==1)? TPoints1:TPoints2;
732 IntPolyh_ArrayOfTriangles &TTriangles=(SurfID==1)? TTriangles1:TTriangles2;
733 Standard_Integer NbSamplesU=(SurfID==1)? NbSamplesU1:NbSamplesU2;
734 Standard_Integer NbSamplesV=(SurfID==1)? NbSamplesV1:NbSamplesV2;
736 TTriangles.Init(2*(NbSamplesU-1)*(NbSamplesV-1));
737 //To provide recursion, I associate a point with two triangles
738 for(Standard_Integer BoucleMeshU=0; BoucleMeshU<NbSamplesU-1; BoucleMeshU++){
739 for(Standard_Integer BoucleMeshV=0; BoucleMeshV<NbSamplesV-1;BoucleMeshV++){
742 TTriangles[CpteurTabTriangles].SetFirstPoint(PntInit); // U V
743 TTriangles[CpteurTabTriangles].SetSecondPoint(PntInit+1); // U V+1
744 TTriangles[CpteurTabTriangles].SetThirdPoint(PntInit+NbSamplesV+1); // U+1 V+1
746 // IF ITS EDGE CONTACTS WITH THE COMMON BOX IP REMAINS = A 1
747 if( ( (TPoints[PntInit].PartOfCommon()) & (TPoints[PntInit+1].PartOfCommon()) )
748 &&( (TPoints[PntInit+1].PartOfCommon()) & (TPoints[PntInit+NbSamplesV+1].PartOfCommon()))
749 &&( (TPoints[PntInit+NbSamplesV+1].PartOfCommon()) & (TPoints[PntInit].PartOfCommon())) )
751 TTriangles[CpteurTabTriangles].SetIndiceIntersectionPossible(0);
753 CpteurTabTriangles++;
756 TTriangles[CpteurTabTriangles].SetFirstPoint(PntInit); // U V
757 TTriangles[CpteurTabTriangles].SetSecondPoint(PntInit+NbSamplesV+1); // U+1 V+1
758 TTriangles[CpteurTabTriangles].SetThirdPoint(PntInit+NbSamplesV); // U+1 V
761 if( ( (TPoints[PntInit].PartOfCommon()) & (TPoints[PntInit+NbSamplesV+1].PartOfCommon()) )
762 &&( (TPoints[PntInit+NbSamplesV+1].PartOfCommon()) & (TPoints[PntInit+NbSamplesV].PartOfCommon()))
763 &&( (TPoints[PntInit+NbSamplesV].PartOfCommon()) & (TPoints[PntInit].PartOfCommon())) )
764 TTriangles[CpteurTabTriangles].SetIndiceIntersectionPossible(0);
767 CpteurTabTriangles++;
769 PntInit++;//Pass to the next point
771 PntInit++;//Pass the last point of the column
773 TTriangles.SetNbItems(CpteurTabTriangles);
774 const Standard_Integer FinTT = TTriangles.NbItems();
778 //=======================================================================
779 //function : LinkEdges2Triangles
780 //purpose : fill the edge fields in Triangle object for the
781 // two array of triangles.
782 //=======================================================================
783 void IntPolyh_MaillageAffinage::LinkEdges2Triangles()
785 const Standard_Integer FinTT1 = TTriangles1.NbItems();
786 const Standard_Integer FinTT2 = TTriangles2.NbItems();
788 for(Standard_Integer uiui1=0; uiui1<FinTT1; uiui1++) {
789 IntPolyh_Triangle & MyTriangle1=TTriangles1[uiui1];
790 if ( (MyTriangle1.FirstEdge()) == -1 ) {
791 MyTriangle1.SetEdgeandOrientation(1,TEdges1);
792 MyTriangle1.SetEdgeandOrientation(2,TEdges1);
793 MyTriangle1.SetEdgeandOrientation(3,TEdges1);
796 for(Standard_Integer uiui2=0; uiui2<FinTT2; uiui2++) {
797 IntPolyh_Triangle & MyTriangle2=TTriangles2[uiui2];
798 if ( (MyTriangle2.FirstEdge()) == -1 ) {
799 MyTriangle2.SetEdgeandOrientation(1,TEdges2);
800 MyTriangle2.SetEdgeandOrientation(2,TEdges2);
801 MyTriangle2.SetEdgeandOrientation(3,TEdges2);
805 //=======================================================================
806 //function : CommonPartRefinement
807 //purpose : Refine systematicaly all marked triangles of both surfaces
808 // REFINING OF THE COMMON
809 //=======================================================================
810 void IntPolyh_MaillageAffinage::CommonPartRefinement()
812 Standard_Integer FinInit1 = TTriangles1.NbItems();
813 for(Standard_Integer i=0; i<FinInit1; i++) {
814 if(TTriangles1[i].IndiceIntersectionPossible()!=0)
815 TTriangles1[i].MiddleRefinement(i,MaSurface1,TPoints1,TTriangles1,TEdges1);
818 Standard_Integer FinInit2=TTriangles2.NbItems();
819 for(Standard_Integer ii=0; ii<FinInit2; ii++) {
820 if(TTriangles2[ii].IndiceIntersectionPossible()!=0)
821 TTriangles2[ii].MiddleRefinement(ii,MaSurface2,TPoints2,TTriangles2,TEdges2);
825 //=======================================================================
826 //function : LocalSurfaceRefinement
827 //purpose : Refine systematicaly all marked triangles of ONE surface
828 //=======================================================================
829 void IntPolyh_MaillageAffinage::LocalSurfaceRefinement(const Standard_Integer SurfID) {
830 //refine locally, but systematically the chosen surface
832 const Standard_Integer FinInit1 = TTriangles1.NbItems();
833 for(Standard_Integer i=0; i<FinInit1; i++) {
834 if(TTriangles1[i].IndiceIntersectionPossible()!=0)
835 TTriangles1[i].MiddleRefinement(i,MaSurface1,TPoints1,TTriangles1,TEdges1);
840 const Standard_Integer FinInit2 = TTriangles2.NbItems();
841 for(Standard_Integer ii=0; ii<FinInit2; ii++) {
842 if(TTriangles2[ii].IndiceIntersectionPossible()!=0)
843 TTriangles2[ii].MiddleRefinement(ii,MaSurface2,TPoints2,TTriangles2,TEdges2);
847 //=======================================================================
848 //function : ComputeDeflections
849 //purpose : Compute deflection for all triangles of one
850 // surface,and sort min and max of deflections
852 // Calculation of the deflection of all triangles
853 // --> deflection max
854 // --> deflection min
855 //=======================================================================
856 void IntPolyh_MaillageAffinage::ComputeDeflections
857 (const Standard_Integer SurfID)
859 Handle(Adaptor3d_HSurface) MaSurface=(SurfID==1)? MaSurface1:MaSurface2;
860 IntPolyh_ArrayOfPoints &TPoints=(SurfID==1)? TPoints1:TPoints2;
861 IntPolyh_ArrayOfTriangles &TTriangles=(SurfID==1)? TTriangles1:TTriangles2;
862 Standard_Real &FlecheMin=(SurfID==1)? FlecheMin1:FlecheMin2;
863 Standard_Real &FlecheMoy=(SurfID==1)? FlecheMoy1:FlecheMoy2;
864 Standard_Real &FlecheMax=(SurfID==1)? FlecheMax1:FlecheMax2;
866 Standard_Integer CpteurTabFleche=0;
867 FlecheMax=-RealLast();
868 FlecheMin=RealLast();
870 const Standard_Integer FinTT = TTriangles.NbItems();
872 for(CpteurTabFleche=0; CpteurTabFleche<FinTT; CpteurTabFleche++) {
873 IntPolyh_Triangle &Triangle = TTriangles[CpteurTabFleche];
874 if ( Triangle.GetFleche() < 0) { //pas normal
878 Triangle.TriangleDeflection(MaSurface, TPoints);
879 Standard_Real Fleche=Triangle.GetFleche();
881 if (Fleche > FlecheMax)
883 if (Fleche < FlecheMin)
888 //=======================================================================
889 //function : TrianglesDeflectionsRefinementBSB
890 //purpose : Refine both surfaces using BoundSortBox as --
891 // rejection. The criterions used to refine a --
892 // triangle are: The deflection The size of the --
893 // bounding boxes (one surface may be very small
894 // compared to the other)
895 //=======================================================================
896 void IntPolyh_MaillageAffinage::TrianglesDeflectionsRefinementBSB()
898 const Standard_Integer FinTT1 = TTriangles1.NbItems();
899 const Standard_Integer FinTT2 = TTriangles2.NbItems();
901 ComputeDeflections(1);
902 // To estimate a surface in general it can be interesting
903 //to calculate all deflections
904 //-- Check deflection at output
906 Standard_Real FlecheCritique1;
907 if(FlecheMin1>FlecheMax1) {
910 else {//fleche min + (flechemax-flechemin) * 80/100
911 FlecheCritique1 = FlecheMin1*0.2+FlecheMax1*0.8;
914 ComputeDeflections(2);
915 //-- Check arrows at output
917 Standard_Real FlecheCritique2;
918 if(FlecheMin2>FlecheMax2) {
922 else {//fleche min + (flechemax-flechemin) * 80/100
923 FlecheCritique2 = FlecheMin2*0.2+FlecheMax2*0.8;
927 Bnd_BoundSortBox BndBSB;
928 Standard_Real diag1,diag2;
929 Standard_Real x0,y0,z0,x1,y1,z1;
931 //The greatest of two bounding boxes created in FillArrayOfPoints is found.
932 //Then this value is weighted depending on the discretization
933 //(NbSamplesU and NbSamplesV)
934 MyBox1.Get(x0,y0,z0,x1,y1,z1);
935 x0-=x1; y0-=y1; z0-=z1;
936 diag1=x0*x0+y0*y0+z0*z0;
937 const Standard_Real NbSamplesUV1=Standard_Real(NbSamplesU1) * Standard_Real(NbSamplesV1);
940 MyBox2.Get(x0,y0,z0,x1,y1,z1);
941 x0-=x1; y0-=y1; z0-=z1;
942 diag2=x0*x0+y0*y0+z0*z0;
943 const Standard_Real NbSamplesUV2=Standard_Real(NbSamplesU2) * Standard_Real(NbSamplesV2);
946 //-- The surface with the greatest bounding box is "discretized"
948 //Standard_Integer NbInterTentees=0;
952 if(FlecheCritique2<diag1) {//the corresponding sizes are not too disproportional
954 Handle(Bnd_HArray1OfBox) HBnd = new Bnd_HArray1OfBox(1,FinTT2);
956 for(Standard_Integer i=0; i<FinTT2; i++){
957 if (TTriangles2[i].IndiceIntersectionPossible()!=0) {
959 const IntPolyh_Triangle& T=TTriangles2[i];
960 const IntPolyh_Point& PA=TPoints2[T.FirstPoint()];
961 const IntPolyh_Point& PB=TPoints2[T.SecondPoint()];
962 const IntPolyh_Point& PC=TPoints2[T.ThirdPoint()];
963 gp_Pnt pntA(PA.X(),PA.Y(),PA.Z());
964 gp_Pnt pntB(PB.X(),PB.Y(),PB.Z());
965 gp_Pnt pntC(PC.X(),PC.Y(),PC.Z());
966 b.Add(pntA);//Box b, which contains triangle i of surface 2 is created./
969 b.Enlarge(T.GetFleche()+MyTolerance);
970 HBnd->SetValue(i+1,b);//Box b is added in the array HBnd
974 //Inititalization of the boundary, sorting of boxes
975 BndBSB.Initialize(HBnd);//contains boxes of 2
977 Standard_Integer FinTT1Init=FinTT1;
978 for(Standard_Integer i_S1=0; i_S1<FinTT1Init; i_S1++) {
979 if(TTriangles1[i_S1].IndiceIntersectionPossible()!=0) {
980 //-- Loop on the boxes of mesh 1
982 const IntPolyh_Triangle& T=TTriangles1[i_S1];
983 const IntPolyh_Point& PA=TPoints1[T.FirstPoint()];
984 const IntPolyh_Point& PB=TPoints1[T.SecondPoint()];
985 const IntPolyh_Point& PC=TPoints1[T.ThirdPoint()];
986 gp_Pnt pntA(PA.X(),PA.Y(),PA.Z());
987 gp_Pnt pntB(PB.X(),PB.Y(),PB.Z());
988 gp_Pnt pntC(PC.X(),PC.Y(),PC.Z());
992 b.Enlarge(T.GetFleche());
993 //-- List of boxes of 2, which touch this box (of 1)
994 const TColStd_ListOfInteger& ListeOf2 = BndBSB.Compare(b);
996 if((ListeOf2.IsEmpty())==0) {
997 IntPolyh_Triangle &Triangle1 = TTriangles1[i_S1];
998 if(Triangle1.GetFleche()>FlecheCritique1)
999 Triangle1.MiddleRefinement(i_S1,MaSurface1,TPoints1,
1000 TTriangles1, TEdges1);
1002 for (TColStd_ListIteratorOfListOfInteger Iter(ListeOf2);
1005 Standard_Integer i_S2=Iter.Value()-1;
1006 //if the box of s1 contacts with the boxes of s2
1007 //the arrow of the triangle is checked
1008 IntPolyh_Triangle & Triangle2 = TTriangles2[i_S2];
1009 if(Triangle2.IndiceIntersectionPossible()!=0)
1010 if(Triangle2.GetFleche()>FlecheCritique2)
1011 Triangle2.MiddleRefinement( i_S2, MaSurface2, TPoints2,
1012 TTriangles2, TEdges2);
1019 //--------------------------------------------------------------------
1020 //FlecheCritique2 > diag1
1024 Handle(Bnd_HArray1OfBox) HBnd = new Bnd_HArray1OfBox(1,FinTT2);
1026 for(Standard_Integer i=0; i<FinTT2; i++){
1027 if (TTriangles2[i].IndiceIntersectionPossible()!=0) {
1029 const IntPolyh_Triangle& T=TTriangles2[i];
1030 const IntPolyh_Point& PA=TPoints2[T.FirstPoint()];
1031 const IntPolyh_Point& PB=TPoints2[T.SecondPoint()];
1032 const IntPolyh_Point& PC=TPoints2[T.ThirdPoint()];
1033 gp_Pnt pntA(PA.X(),PA.Y(),PA.Z());
1034 gp_Pnt pntB(PB.X(),PB.Y(),PB.Z());
1035 gp_Pnt pntC(PC.X(),PC.Y(),PC.Z());
1036 b.Add(pntA);//Box b, which contains triangle i of surface 2 is created/
1039 b.Enlarge(T.GetFleche()+MyTolerance);
1040 //-- BndBSB.Add(b,i+1);
1041 HBnd->SetValue(i+1,b);//Box b is added in array HBnd
1045 //Inititalization of the ouput bounding box
1046 BndBSB.Initialize(HBnd);//contains boxes of 2
1049 //The bounding box Be1 of surface1 is compared BSB of surface2
1050 const TColStd_ListOfInteger& ListeOf2 = BndBSB.Compare(MyBox1);
1052 if((ListeOf2.IsEmpty())==0) {
1053 //if the bounding box Be1 of s1 contacts with
1054 //the boxes of s2 the deflection of triangle of s2 is checked
1056 // Be1 is very small in relation to Be2
1057 //The criterion of refining for surface2 depends on the size of Be1
1058 //As it is known that this criterion should be minimized,
1059 //the smallest side of the bounding box is taken
1060 MyBox1.Get(x0,y0,z0,x1,y1,z1);
1061 Standard_Real dx=Abs(x1-x0);
1062 Standard_Real dy=Abs(y1-y0);
1063 Standard_Real diag=Abs(z1-z0);
1064 Standard_Real dd=-1.0;
1069 if (diag>dd) diag=dd;
1071 //if Be1 contacts with the boxes of s2, the deflection
1072 //of the triangles of s2 is checked (greater)
1073 //in relation to the size of Be1 (smaller)
1074 for (TColStd_ListIteratorOfListOfInteger Iter(ListeOf2);
1077 Standard_Integer i_S2=Iter.Value()-1;
1079 IntPolyh_Triangle & Triangle2=TTriangles2[i_S2];
1080 if(Triangle2.IndiceIntersectionPossible()) {
1082 //calculation of the criterion of refining
1083 //The deflection of the greater is compared to the size of the smaller
1084 Standard_Real CritereAffinage=0.0;
1085 Standard_Real DiagPonderation=0.5;
1086 CritereAffinage = diag*DiagPonderation;
1087 if(Triangle2.GetFleche()>CritereAffinage)
1088 Triangle2.MultipleMiddleRefinement2(CritereAffinage, MyBox1, i_S2,
1089 MaSurface2, TPoints2,
1090 TTriangles2,TEdges2);
1092 else Triangle2.MiddleRefinement(i_S2,MaSurface2,TPoints2,
1093 TTriangles2, TEdges2);
1101 else { //-- The greater is discretised
1103 if(FlecheCritique1<diag2) {//the respective sizes are not to much disproportional
1105 Handle(Bnd_HArray1OfBox) HBnd = new Bnd_HArray1OfBox(1,FinTT1);
1107 for(Standard_Integer i=0; i<FinTT1; i++){
1108 if(TTriangles1[i].IndiceIntersectionPossible()!=0) {
1110 const IntPolyh_Triangle& T=TTriangles1[i];
1111 const IntPolyh_Point& PA=TPoints1[T.FirstPoint()];
1112 const IntPolyh_Point& PB=TPoints1[T.SecondPoint()];
1113 const IntPolyh_Point& PC=TPoints1[T.ThirdPoint()];
1114 gp_Pnt pntA(PA.X(),PA.Y(),PA.Z());
1115 gp_Pnt pntB(PB.X(),PB.Y(),PB.Z());
1116 gp_Pnt pntC(PC.X(),PC.Y(),PC.Z());
1117 b.Add(pntA);//Box b, which contains triangle i of surface 2 is created.
1120 b.Enlarge(T.GetFleche()+MyTolerance);
1121 HBnd->SetValue(i+1,b);//Boite b is added in the array HBnd
1124 BndBSB.Initialize(HBnd);
1126 Standard_Integer FinTT2init=FinTT2;
1127 for(Standard_Integer i_S2=0; i_S2<FinTT2init; i_S2++) {
1128 if (TTriangles2[i_S2].IndiceIntersectionPossible()!=0) {
1129 //-- Loop on the boxes of mesh 2
1131 const IntPolyh_Triangle& T=TTriangles2[i_S2];
1132 const IntPolyh_Point& PA=TPoints2[T.FirstPoint()];
1133 const IntPolyh_Point& PB=TPoints2[T.SecondPoint()];
1134 const IntPolyh_Point& PC=TPoints2[T.ThirdPoint()];
1135 gp_Pnt pntA(PA.X(),PA.Y(),PA.Z());
1136 gp_Pnt pntB(PB.X(),PB.Y(),PB.Z());
1137 gp_Pnt pntC(PC.X(),PC.Y(),PC.Z());
1141 b.Enlarge(T.GetFleche()+MyTolerance);
1142 //-- List of boxes of 1 touching this box (of 2)
1143 const TColStd_ListOfInteger& ListeOf1 = BndBSB.Compare(b);
1144 IntPolyh_Triangle & Triangle2=TTriangles2[i_S2];
1145 if((ListeOf1.IsEmpty())==0) {
1147 if(Triangle2.GetFleche()>FlecheCritique2)
1148 Triangle2.MiddleRefinement(i_S2,MaSurface2,TPoints2,
1149 TTriangles2, TEdges2);
1151 for (TColStd_ListIteratorOfListOfInteger Iter(ListeOf1);
1154 Standard_Integer i_S1=Iter.Value()-1;
1155 IntPolyh_Triangle & Triangle1=TTriangles1[i_S1];
1156 if (Triangle1.IndiceIntersectionPossible())
1157 if(Triangle1.GetFleche()>FlecheCritique1)
1158 Triangle1.MiddleRefinement(i_S1,MaSurface1,TPoints1,
1159 TTriangles1, TEdges1);
1165 //-----------------------------------------------------------------------------
1166 else {// FlecheCritique1>diag2
1169 Handle(Bnd_HArray1OfBox) HBnd = new Bnd_HArray1OfBox(1,FinTT1);
1171 for(Standard_Integer i=0; i<FinTT1; i++){
1172 if (TTriangles1[i].IndiceIntersectionPossible()!=0) {
1174 const IntPolyh_Triangle& T=TTriangles1[i];
1175 const IntPolyh_Point& PA=TPoints1[T.FirstPoint()];
1176 const IntPolyh_Point& PB=TPoints1[T.SecondPoint()];
1177 const IntPolyh_Point& PC=TPoints1[T.ThirdPoint()];
1178 gp_Pnt pntA(PA.X(),PA.Y(),PA.Z());
1179 gp_Pnt pntB(PB.X(),PB.Y(),PB.Z());
1180 gp_Pnt pntC(PC.X(),PC.Y(),PC.Z());
1181 b.Add(pntA);//Box b, which contains triangle i of surface 1 is created./
1184 b.Enlarge(T.GetFleche()+MyTolerance);
1185 HBnd->SetValue(i+1,b);//Box b is added in the array HBnd
1189 //Inititalisation of the boundary output box
1190 BndBSB.Initialize(HBnd);//contains boxes of 1
1192 //Bounding box Be2 of surface2 is compared to BSB of surface1
1193 const TColStd_ListOfInteger& ListeOf1 = BndBSB.Compare(MyBox2);
1195 if((ListeOf1.IsEmpty())==0) {
1196 //if the bounding box Be2 of s2 contacts
1197 //with boxes of s1 the deflection of the triangle of s1 is checked
1199 // Be2 is very small compared to Be1
1200 //The criterion of refining for surface1 depends on the size of Be2
1201 //As this criterion should be minimized,
1202 //the smallest side of the bounding box is taken
1203 MyBox2.Get(x0,y0,z0,x1,y1,z1);
1204 Standard_Real dx=Abs(x1-x0);
1205 Standard_Real dy=Abs(y1-y0);
1206 Standard_Real diag=Abs(z1-z0);
1207 Standard_Real dd=-1.0;
1212 if (diag>dd) diag=dd;
1214 //if Be2 contacts with boxes of s1, the deflection of
1215 //triangles of s1 (greater) is checked
1216 //comparatively to the size of Be2 (smaller).
1217 for (TColStd_ListIteratorOfListOfInteger Iter(ListeOf1);
1220 Standard_Integer i_S1=Iter.Value()-1;
1222 IntPolyh_Triangle & Triangle1=TTriangles1[i_S1];
1223 if(Triangle1.IndiceIntersectionPossible()) {
1225 //calculation of the criterion of refining
1226 //The deflection of the greater is compared
1227 //with the size of the smaller.
1228 Standard_Real CritereAffinage=0.0;
1229 Standard_Real DiagPonderation=0.5;
1230 CritereAffinage = diag*DiagPonderation;;
1231 if(Triangle1.GetFleche()>CritereAffinage)
1232 Triangle1.MultipleMiddleRefinement2(CritereAffinage,MyBox2, i_S1,
1233 MaSurface1, TPoints1,
1234 TTriangles1, TEdges1);
1236 else Triangle1.MiddleRefinement(i_S1,MaSurface1,TPoints1,
1237 TTriangles1, TEdges1);
1245 //=======================================================================
1247 //purpose : This function is used for the function project6
1248 //=======================================================================
1249 inline Standard_Real maxSR(const Standard_Real a,
1250 const Standard_Real b,
1251 const Standard_Real c)
1253 Standard_Real t = a;
1258 //=======================================================================
1260 //purpose : This function is used for the function project6
1261 //=======================================================================
1262 inline Standard_Real minSR(const Standard_Real a,
1263 const Standard_Real b,
1264 const Standard_Real c)
1266 Standard_Real t = a;
1271 //=======================================================================
1272 //function : project6
1273 //purpose : This function is used for the function TriContact
1274 //=======================================================================
1275 Standard_Integer project6(const IntPolyh_Point &ax,
1276 const IntPolyh_Point &p1,
1277 const IntPolyh_Point &p2,
1278 const IntPolyh_Point &p3,
1279 const IntPolyh_Point &q1,
1280 const IntPolyh_Point &q2,
1281 const IntPolyh_Point &q3)
1283 Standard_Real P1 = ax.Dot(p1);
1284 Standard_Real P2 = ax.Dot(p2);
1285 Standard_Real P3 = ax.Dot(p3);
1286 Standard_Real Q1 = ax.Dot(q1);
1287 Standard_Real Q2 = ax.Dot(q2);
1288 Standard_Real Q3 = ax.Dot(q3);
1290 Standard_Real mx1 = maxSR(P1, P2, P3);
1291 Standard_Real mn1 = minSR(P1, P2, P3);
1292 Standard_Real mx2 = maxSR(Q1, Q2, Q3);
1293 Standard_Real mn2 = minSR(Q1, Q2, Q3);
1295 if (mn1 > mx2) return 0;
1296 if (mn2 > mx1) return 0;
1299 //=======================================================================
1300 //function : TriContact
1301 //purpose : This fonction Check if two triangles are in
1302 // contact or no, return 1 if yes, return 0
1304 //=======================================================================
1305 Standard_Integer IntPolyh_MaillageAffinage::TriContact
1306 (const IntPolyh_Point &P1,
1307 const IntPolyh_Point &P2,
1308 const IntPolyh_Point &P3,
1309 const IntPolyh_Point &Q1,
1310 const IntPolyh_Point &Q2,
1311 const IntPolyh_Point &Q3,
1312 Standard_Real &Angle) const
1315 The first triangle is (p1,p2,p3). The other is (q1,q2,q3).
1316 The edges are (e1,e2,e3) and (f1,f2,f3).
1317 The normals are n1 and m1
1318 Outwards are (g1,g2,g3) and (h1,h2,h3).*/
1320 if(maxSR(P1.X(),P2.X(),P3.X())<minSR(Q1.X(),Q2.X(),Q3.X())) return(0);
1321 if(maxSR(P1.Y(),P2.Y(),P3.Y())<minSR(Q1.Y(),Q2.Y(),Q3.Y())) return(0);
1322 if(maxSR(P1.Z(),P2.Z(),P3.Z())<minSR(Q1.Z(),Q2.Z(),Q3.Z())) return(0);
1324 if(minSR(P1.X(),P2.X(),P3.X())>maxSR(Q1.X(),Q2.X(),Q3.X())) return(0);
1325 if(minSR(P1.Y(),P2.Y(),P3.Y())>maxSR(Q1.Y(),Q2.Y(),Q3.Y())) return(0);
1326 if(minSR(P1.Z(),P2.Z(),P3.Z())>maxSR(Q1.Z(),Q2.Z(),Q3.Z())) return(0);
1328 IntPolyh_Point p1, p2, p3;
1329 IntPolyh_Point q1, q2, q3;
1330 IntPolyh_Point e1, e2, e3;
1331 IntPolyh_Point f1, f2, f3;
1332 IntPolyh_Point g1, g2, g3;
1333 IntPolyh_Point h1, h2, h3;
1334 IntPolyh_Point n1, m1;
1337 IntPolyh_Point ef11, ef12, ef13;
1338 IntPolyh_Point ef21, ef22, ef23;
1339 IntPolyh_Point ef31, ef32, ef33;
1341 z.SetX(0.0); z.SetY(0.0); z.SetZ(0.0);
1344 p1.SetX(P1.X() - P1.X()); p1.SetY(P1.Y() - P1.Y()); p1.SetZ(P1.Z() - P1.Z());
1345 p2.SetX(P2.X() - P1.X()); p2.SetY(P2.Y() - P1.Y()); p2.SetZ(P2.Z() - P1.Z());
1346 p3.SetX(P3.X() - P1.X()); p3.SetY(P3.Y() - P1.Y()); p3.SetZ(P3.Z() - P1.Z());
1348 q1.SetX(Q1.X() - P1.X()); q1.SetY(Q1.Y() - P1.Y()); q1.SetZ(Q1.Z() - P1.Z());
1349 q2.SetX(Q2.X() - P1.X()); q2.SetY(Q2.Y() - P1.Y()); q2.SetZ(Q2.Z() - P1.Z());
1350 q3.SetX(Q3.X() - P1.X()); q3.SetY(Q3.Y() - P1.Y()); q3.SetZ(Q3.Z() - P1.Z());
1352 e1.SetX(p2.X() - p1.X()); e1.SetY(p2.Y() - p1.Y()); e1.SetZ(p2.Z() - p1.Z());
1353 e2.SetX(p3.X() - p2.X()); e2.SetY(p3.Y() - p2.Y()); e2.SetZ(p3.Z() - p2.Z());
1354 e3.SetX(p1.X() - p3.X()); e3.SetY(p1.Y() - p3.Y()); e3.SetZ(p1.Z() - p3.Z());
1356 f1.SetX(q2.X() - q1.X()); f1.SetY(q2.Y() - q1.Y()); f1.SetZ(q2.Z() - q1.Z());
1357 f2.SetX(q3.X() - q2.X()); f2.SetY(q3.Y() - q2.Y()); f2.SetZ(q3.Z() - q2.Z());
1358 f3.SetX(q1.X() - q3.X()); f3.SetY(q1.Y() - q3.Y()); f3.SetZ(q1.Z() - q3.Z());
1360 n1.Cross(e1, e2); //normal to the first triangle
1361 m1.Cross(f1, f2); //normal to the second triangle
1380 // Now the testing is done
1382 if (!project6(n1, p1, p2, p3, q1, q2, q3)) return 0; //T2 is not higher or lower than T1
1383 if (!project6(m1, p1, p2, p3, q1, q2, q3)) return 0; //T1 is not higher of lower than T2
1385 if (!project6(ef11, p1, p2, p3, q1, q2, q3)) return 0;
1386 if (!project6(ef12, p1, p2, p3, q1, q2, q3)) return 0;
1387 if (!project6(ef13, p1, p2, p3, q1, q2, q3)) return 0;
1388 if (!project6(ef21, p1, p2, p3, q1, q2, q3)) return 0;
1389 if (!project6(ef22, p1, p2, p3, q1, q2, q3)) return 0;
1390 if (!project6(ef23, p1, p2, p3, q1, q2, q3)) return 0;
1391 if (!project6(ef31, p1, p2, p3, q1, q2, q3)) return 0;
1392 if (!project6(ef32, p1, p2, p3, q1, q2, q3)) return 0;
1393 if (!project6(ef33, p1, p2, p3, q1, q2, q3)) return 0;
1395 if (!project6(g1, p1, p2, p3, q1, q2, q3)) return 0; //T2 is outside of T1 in the plane of T1
1396 if (!project6(g2, p1, p2, p3, q1, q2, q3)) return 0; //T2 is outside of T1 in the plane of T1
1397 if (!project6(g3, p1, p2, p3, q1, q2, q3)) return 0; //T2 is outside of T1 in the plane of T1
1398 if (!project6(h1, p1, p2, p3, q1, q2, q3)) return 0; //T1 is outside of T2 in the plane of T2
1399 if (!project6(h2, p1, p2, p3, q1, q2, q3)) return 0; //T1 is outside of T2 in the plane of T2
1400 if (!project6(h3, p1, p2, p3, q1, q2, q3)) return 0; //T1 is outside of T2 in the plane of T2
1402 //Calculation of cosinus angle between two normals
1403 Standard_Real SqModn1=-1.0;
1404 Standard_Real SqModm1=-1.0;
1405 SqModn1=n1.SquareModulus();
1406 if (SqModn1>SquareMyConfusionPrecision){
1407 SqModm1=m1.SquareModulus();
1409 if (SqModm1>SquareMyConfusionPrecision) {
1410 Angle=(n1.Dot(m1))/(sqrt(SqModn1)*sqrt(SqModm1));
1414 //=======================================================================
1415 //function : TestNbPoints
1416 //purpose : This function is used by StartingPointsResearch() to control
1417 // the number of points found keep the result in conformity (1 or 2 points)
1418 // void TestNbPoints(const Standard_Integer TriSurfID,
1419 //=======================================================================
1420 void TestNbPoints(const Standard_Integer ,
1421 Standard_Integer &NbPoints,
1422 Standard_Integer &NbPointsTotal,
1423 const IntPolyh_StartPoint &Pt1,
1424 const IntPolyh_StartPoint &Pt2,
1425 IntPolyh_StartPoint &SP1,
1426 IntPolyh_StartPoint &SP2)
1428 // already checked in TriangleEdgeContact2
1429 // if( (NbPoints==2)&&(Pt1.CheckSameSP(Pt2)) ) NbPoints=1;
1435 if ( (NbPoints==1)&&(NbPointsTotal==0) ) {
1439 else if ( (NbPoints==1)&&(NbPointsTotal==1) ) {
1440 if(Pt1.CheckSameSP(SP1)!=1) {
1445 else if( (NbPoints==1)&&(NbPointsTotal==2) ) {
1446 if ( (SP1.CheckSameSP(Pt1))||(SP2.CheckSameSP(Pt1)) )
1448 else NbPointsTotal=3;
1450 else if( (NbPoints==2)&&(NbPointsTotal==0) ) {
1455 else if( (NbPoints==2)&&(NbPointsTotal==1) ) {//there is also Pt1 != Pt2
1456 if(SP1.CheckSameSP(Pt1)) {
1460 else if (SP1.CheckSameSP(Pt2)) {
1464 else NbPointsTotal=3;///case SP1!=Pt1 && SP1!=Pt2!
1466 else if( (NbPoints==2)&&(NbPointsTotal==2) ) {//there is also SP1!=SP2
1467 if( (SP1.CheckSameSP(Pt1))||(SP1.CheckSameSP(Pt2)) ) {
1468 if( (SP2.CheckSameSP(Pt1))||(SP2.CheckSameSP(Pt2)) )
1470 else NbPointsTotal=3;
1472 else NbPointsTotal=3;
1476 //=======================================================================
1477 //function : StartingPointsResearch
1479 //=======================================================================
1480 Standard_Integer IntPolyh_MaillageAffinage::StartingPointsResearch
1481 (const Standard_Integer T1,
1482 const Standard_Integer T2,
1483 IntPolyh_StartPoint &SP1,
1484 IntPolyh_StartPoint &SP2) const
1486 const IntPolyh_Point &P1=TPoints1[TTriangles1[T1].FirstPoint()];
1487 const IntPolyh_Point &P2=TPoints1[TTriangles1[T1].SecondPoint()];
1488 const IntPolyh_Point &P3=TPoints1[TTriangles1[T1].ThirdPoint()];
1489 const IntPolyh_Point &Q1=TPoints2[TTriangles2[T2].FirstPoint()];
1490 const IntPolyh_Point &Q2=TPoints2[TTriangles2[T2].SecondPoint()];
1491 const IntPolyh_Point &Q3=TPoints2[TTriangles2[T2].ThirdPoint()];
1494 /* The first triangle is (p1,p2,p3). The other is (q1,q2,q3).
1495 The sides are (e1,e2,e3) and (f1,f2,f3).
1496 The normals are n1 and m1*/
1498 const IntPolyh_Point e1=P2-P1;
1499 const IntPolyh_Point e2=P3-P2;
1500 const IntPolyh_Point e3=P1-P3;
1502 const IntPolyh_Point f1=Q2-Q1;
1503 const IntPolyh_Point f2=Q3-Q2;
1504 const IntPolyh_Point f3=Q1-Q3;
1507 IntPolyh_Point nn1,mm1;
1508 nn1.Cross(e1, e2); //normal of the first triangle
1509 mm1.Cross(f1, f2); //normal of the second triangle
1511 Standard_Real nn1modulus, mm1modulus;
1512 nn1modulus=sqrt(nn1.SquareModulus());
1513 mm1modulus=sqrt(mm1.SquareModulus());
1515 //-------------------------------------------------------
1516 ///calculation of intersection points between two triangles
1517 //-------------------------------------------------------
1518 Standard_Integer NbPoints=0;
1519 Standard_Integer NbPointsTotal=0;
1520 IntPolyh_StartPoint Pt1,Pt2;
1524 if(Abs(nn1modulus)<MyConfusionPrecision){//10.0e-20) {
1528 const IntPolyh_Point n1=nn1.Divide(nn1modulus);
1530 if(NbPointsTotal<2) {
1531 NbPoints=TriangleEdgeContact(1,1,P1,P2,P3,e1,e2,e3,Q1,Q2,f1,n1,Pt1,Pt2);
1532 TestNbPoints(1,NbPoints,NbPointsTotal,Pt1,Pt2,SP1,SP2);
1535 if(NbPointsTotal<2) {
1536 NbPoints=TriangleEdgeContact(1,2,P1,P2,P3,e1,e2,e3,Q2,Q3,f2,n1,Pt1,Pt2);
1537 TestNbPoints(1,NbPoints,NbPointsTotal,Pt1,Pt2,SP1,SP2);
1540 if(NbPointsTotal<2) {
1541 NbPoints=TriangleEdgeContact(1,3,P1,P2,P3,e1,e2,e3,Q3,Q1,f3,n1,Pt1,Pt2);
1542 TestNbPoints(1,NbPoints,NbPointsTotal,Pt1,Pt2,SP1,SP2);
1547 if(Abs(mm1modulus)<MyConfusionPrecision) {//10.0e-20){
1551 const IntPolyh_Point m1=mm1.Divide(mm1modulus);
1553 if(NbPointsTotal<2) {
1554 NbPoints=TriangleEdgeContact(2,1,Q1,Q2,Q3,f1,f2,f3,P1,P2,e1,m1,Pt1,Pt2);
1555 TestNbPoints(2,NbPoints,NbPointsTotal,Pt1,Pt2,SP1,SP2);
1558 if(NbPointsTotal<2) {
1559 NbPoints=TriangleEdgeContact(2,2,Q1,Q2,Q3,f1,f2,f3,P2,P3,e2,m1,Pt1,Pt2);
1560 TestNbPoints(2,NbPoints,NbPointsTotal,Pt1,Pt2,SP1,SP2);
1563 if(NbPointsTotal<2) {
1564 NbPoints=TriangleEdgeContact(2,3,Q1,Q2,Q3,f1,f2,f3,P3,P1,e3,m1,Pt1,Pt2);
1565 TestNbPoints(2,NbPoints,NbPointsTotal,Pt1,Pt2,SP1,SP2);
1569 /* if( (NbPointsTotal >1)&&( Abs(SP1.U1()-SP2.U1())<MyConfusionPrecision)
1570 &&(Abs(SP1.V1()-SP2.V1())<MyConfusionPrecision) )*/
1571 if( (NbPoints)&&(SP1.CheckSameSP(SP2)) )
1574 SP1.SetCoupleValue(T1,T2);
1575 SP2.SetCoupleValue(T1,T2);
1576 return (NbPointsTotal);
1578 //=======================================================================
1579 //function : StartingPointsResearch2
1580 //purpose : From two triangles compute intersection points.
1581 // If I found more than two intersection points
1582 // it means that those triangle are coplanar
1583 //=======================================================================
1584 Standard_Integer IntPolyh_MaillageAffinage::StartingPointsResearch2
1585 (const Standard_Integer T1,
1586 const Standard_Integer T2,
1587 IntPolyh_StartPoint &SP1,
1588 IntPolyh_StartPoint &SP2) const
1590 const IntPolyh_Triangle &Tri1=TTriangles1[T1];
1591 const IntPolyh_Triangle &Tri2=TTriangles2[T2];
1593 const IntPolyh_Point &P1=TPoints1[Tri1.FirstPoint()];
1594 const IntPolyh_Point &P2=TPoints1[Tri1.SecondPoint()];
1595 const IntPolyh_Point &P3=TPoints1[Tri1.ThirdPoint()];
1596 const IntPolyh_Point &Q1=TPoints2[Tri2.FirstPoint()];
1597 const IntPolyh_Point &Q2=TPoints2[Tri2.SecondPoint()];
1598 const IntPolyh_Point &Q3=TPoints2[Tri2.ThirdPoint()];
1602 /* The first triangle is (p1,p2,p3). The other is (q1,q2,q3).
1603 The sides are (e1,e2,e3) and (f1,f2,f3).
1604 The normals are n1 and m1*/
1606 const IntPolyh_Point e1=P2-P1;
1607 const IntPolyh_Point e2=P3-P2;
1608 const IntPolyh_Point e3=P1-P3;
1610 const IntPolyh_Point f1=Q2-Q1;
1611 const IntPolyh_Point f2=Q3-Q2;
1612 const IntPolyh_Point f3=Q1-Q3;
1615 IntPolyh_Point nn1,mm1;
1616 nn1.Cross(e1, e2); //normal to the first triangle
1617 mm1.Cross(f1, f2); //normal to the second triangle
1619 Standard_Real nn1modulus, mm1modulus;
1620 nn1modulus=sqrt(nn1.SquareModulus());
1621 mm1modulus=sqrt(mm1.SquareModulus());
1623 //-------------------------------------------------
1624 ///calculation of intersections points between triangles
1625 //-------------------------------------------------
1626 Standard_Integer NbPoints=0;
1627 Standard_Integer NbPointsTotal=0;
1631 if(Abs(nn1modulus)<MyConfusionPrecision) {//10.0e-20){
1635 const IntPolyh_Point n1=nn1.Divide(nn1modulus);
1637 if(NbPointsTotal<3) {
1638 IntPolyh_StartPoint Pt1,Pt2;
1639 NbPoints=TriangleEdgeContact2(1,1,Tri1,Tri2,P1,P2,P3,e1,e2,e3,Q1,Q2,f1,n1,Pt1,Pt2);
1640 TestNbPoints(1,NbPoints,NbPointsTotal,Pt1,Pt2,SP1,SP2);
1643 if(NbPointsTotal<3) {
1644 IntPolyh_StartPoint Pt1,Pt2;
1645 NbPoints=TriangleEdgeContact2(1,2,Tri1,Tri2,P1,P2,P3,e1,e2,e3,Q2,Q3,f2,n1,Pt1,Pt2);
1646 TestNbPoints(1,NbPoints,NbPointsTotal,Pt1,Pt2,SP1,SP2);
1649 if(NbPointsTotal<3) {
1650 IntPolyh_StartPoint Pt1,Pt2;
1651 NbPoints=TriangleEdgeContact2(1,3,Tri1,Tri2,P1,P2,P3,e1,e2,e3,Q3,Q1,f3,n1,Pt1,Pt2);
1652 TestNbPoints(1,NbPoints,NbPointsTotal,Pt1,Pt2,SP1,SP2);
1657 if(Abs(mm1modulus)<MyConfusionPrecision) {//10.0e-20){
1661 const IntPolyh_Point m1=mm1.Divide(mm1modulus);
1663 if(NbPointsTotal<3) {
1664 IntPolyh_StartPoint Pt1,Pt2;
1665 NbPoints=TriangleEdgeContact2(2,1,Tri1,Tri2,Q1,Q2,Q3,f1,f2,f3,P1,P2,e1,m1,Pt1,Pt2);
1666 TestNbPoints(2,NbPoints,NbPointsTotal,Pt1,Pt2,SP1,SP2);
1669 if(NbPointsTotal<3) {
1670 IntPolyh_StartPoint Pt1,Pt2;
1671 NbPoints=TriangleEdgeContact2(2,2,Tri1,Tri2,Q1,Q2,Q3,f1,f2,f3,P2,P3,e2,m1,Pt1,Pt2);
1672 TestNbPoints(2,NbPoints,NbPointsTotal,Pt1,Pt2,SP1,SP2);
1675 if(NbPointsTotal<3) {
1676 IntPolyh_StartPoint Pt1,Pt2;
1677 NbPoints=TriangleEdgeContact2(2,3,Tri1,Tri2,Q1,Q2,Q3,f1,f2,f3,P3,P1,e3,m1,Pt1,Pt2);
1678 TestNbPoints(2,NbPoints,NbPointsTotal,Pt1,Pt2,SP1,SP2);
1681 // no need already checked in TestNbPoints
1682 // if( (NbPointsTotal==2)&&(SP1.CheckSameSP(SP2)) ) {
1684 //SP1.SetCoupleValue(T1,T2);
1687 if(NbPointsTotal==2) {
1688 SP1.SetCoupleValue(T1,T2);
1689 SP2.SetCoupleValue(T1,T2);
1691 else if(NbPointsTotal==1)
1692 SP1.SetCoupleValue(T1,T2);
1693 else if(NbPointsTotal==3)
1694 SP1.SetCoupleValue(T1,T2);
1696 return (NbPointsTotal);
1698 //=======================================================================
1699 //function : NextStartingPointsResearch
1701 //=======================================================================
1702 Standard_Integer IntPolyh_MaillageAffinage::NextStartingPointsResearch
1703 (const Standard_Integer T1,
1704 const Standard_Integer T2,
1705 const IntPolyh_StartPoint &SPInit,
1706 IntPolyh_StartPoint &SPNext) const
1708 Standard_Integer NbPointsTotal=0;
1709 if( (T1<0)||(T2<0) ) NbPointsTotal=0;
1711 const IntPolyh_Point &P1=TPoints1[TTriangles1[T1].FirstPoint()];
1712 const IntPolyh_Point &P2=TPoints1[TTriangles1[T1].SecondPoint()];
1713 const IntPolyh_Point &P3=TPoints1[TTriangles1[T1].ThirdPoint()];
1714 const IntPolyh_Point &Q1=TPoints2[TTriangles2[T2].FirstPoint()];
1715 const IntPolyh_Point &Q2=TPoints2[TTriangles2[T2].SecondPoint()];
1716 const IntPolyh_Point &Q3=TPoints2[TTriangles2[T2].ThirdPoint()];
1718 /* The first triangle is (p1,p2,p3). The other is (q1,q2,q3).
1719 The sides are (e1,e2,e3) and (f1,f2,f3).
1720 The normals are n1 and m1*/
1722 const IntPolyh_Point e1=P2-P1;
1723 const IntPolyh_Point e2=P3-P2;
1724 const IntPolyh_Point e3=P1-P3;
1726 const IntPolyh_Point f1=Q2-Q1;
1727 const IntPolyh_Point f2=Q3-Q2;
1728 const IntPolyh_Point f3=Q1-Q3;
1730 IntPolyh_Point nn1,mm1;
1731 nn1.Cross(e1, e2); //normal to the first triangle
1732 mm1.Cross(f1, f2); //normal to the second triangle
1734 Standard_Real nn1modulus, mm1modulus;
1735 nn1modulus=sqrt(nn1.SquareModulus());
1736 mm1modulus=sqrt(mm1.SquareModulus());
1738 //-------------------------------------------------
1739 ///calculation of intersections points between triangles
1740 //-------------------------------------------------
1741 Standard_Integer NbPoints=0;
1742 IntPolyh_StartPoint SP1,SP2;
1745 if(Abs(nn1modulus)<MyConfusionPrecision) {//10.0e-20){
1749 const IntPolyh_Point n1=nn1.Divide(nn1modulus);
1751 if(NbPointsTotal<3) {
1752 IntPolyh_StartPoint Pt1,Pt2;
1753 NbPoints=TriangleEdgeContact(1,1,P1,P2,P3,e1,e2,e3,Q1,Q2,f1,n1,Pt1,Pt2);
1754 TestNbPoints(1,NbPoints,NbPointsTotal,Pt1,Pt2,SP1,SP2);
1757 if(NbPointsTotal<3) {
1758 IntPolyh_StartPoint Pt1,Pt2;
1759 NbPoints=TriangleEdgeContact(1,2,P1,P2,P3,e1,e2,e3,Q2,Q3,f2,n1,Pt1,Pt2);
1760 TestNbPoints(1,NbPoints,NbPointsTotal,Pt1,Pt2,SP1,SP2);
1763 if(NbPointsTotal<3) {
1764 IntPolyh_StartPoint Pt1,Pt2;
1765 NbPoints=TriangleEdgeContact(1,3,P1,P2,P3,e1,e2,e3,Q3,Q1,f3,n1,Pt1,Pt2);
1766 TestNbPoints(1,NbPoints,NbPointsTotal,Pt1,Pt2,SP1,SP2);
1771 if(Abs(mm1modulus)<MyConfusionPrecision) {//10.0e-20){
1775 const IntPolyh_Point m1=mm1.Divide(mm1modulus);
1777 if(NbPointsTotal<3) {
1778 IntPolyh_StartPoint Pt1,Pt2;
1779 NbPoints=TriangleEdgeContact(2,1,Q1,Q2,Q3,f1,f2,f3,P1,P2,e1,m1,Pt1,Pt2);
1780 TestNbPoints(2,NbPoints,NbPointsTotal,Pt1,Pt2,SP1,SP2);
1783 if(NbPointsTotal<3) {
1784 IntPolyh_StartPoint Pt1,Pt2;
1785 NbPoints=TriangleEdgeContact(2,2,Q1,Q2,Q3,f1,f2,f3,P2,P3,e2,m1,Pt1,Pt2);
1786 TestNbPoints(2,NbPoints,NbPointsTotal,Pt1,Pt2,SP1,SP2);
1789 if(NbPointsTotal<3) {
1790 IntPolyh_StartPoint Pt1,Pt2;
1791 NbPoints=TriangleEdgeContact(2,3,Q1,Q2,Q3,f1,f2,f3,P3,P1,e3,m1,Pt1,Pt2);
1792 TestNbPoints(2,NbPoints,NbPointsTotal,Pt1,Pt2,SP1,SP2);
1796 if (NbPointsTotal==1) {
1797 /* if( (Abs(SP1.U1()-SPInit.U1())<MyConfusionPrecision)
1798 &&(Abs(SP1.V1()-SPInit.V1())<MyConfusionPrecision) ) */
1799 if(SP1.CheckSameSP(SP2))
1807 // if ( (NbPointsTotal==2)&&( Abs(SP1.U1()-SPInit.U1())<MyConfusionPrecision)
1808 //&&( Abs(SP1.V1()-SPInit.V1())<MyConfusionPrecision) ) {
1809 if( (NbPointsTotal==2)&&(SP1.CheckSameSP(SPInit)) ) {
1810 NbPointsTotal=1;//SP1 et SPInit sont identiques
1813 // if( (NbPointsTotal==2)&&( Abs(SP2.U1()-SPInit.U1())<MyConfusionPrecision)
1814 //&&( Abs(SP2.V1()-SPInit.V1())<MyConfusionPrecision) ) {
1815 if( (NbPointsTotal==2)&&(SP2.CheckSameSP(SPInit)) ) {
1816 NbPointsTotal=1;//SP2 et SPInit sont identiques
1819 if(NbPointsTotal>1) {
1823 SPNext.SetCoupleValue(T1,T2);
1824 return (NbPointsTotal);
1826 //=======================================================================
1827 //function : NextStartingPointsResearch2
1828 //purpose : from two triangles and an intersection point I
1829 // seach the other point (if it exist).
1830 // This function is used by StartPointChain
1831 //=======================================================================
1832 Standard_Integer IntPolyh_MaillageAffinage::NextStartingPointsResearch2
1833 (const Standard_Integer T1,
1834 const Standard_Integer T2,
1835 const IntPolyh_StartPoint &SPInit,
1836 IntPolyh_StartPoint &SPNext) const
1838 Standard_Integer NbPointsTotal=0;
1839 Standard_Integer EdgeInit1=SPInit.E1();
1840 Standard_Integer EdgeInit2=SPInit.E2();
1841 if( (T1<0)||(T2<0) ) NbPointsTotal=0;
1844 const IntPolyh_Triangle &Tri1=TTriangles1[T1];
1845 const IntPolyh_Triangle &Tri2=TTriangles2[T2];
1847 const IntPolyh_Point &P1=TPoints1[Tri1.FirstPoint()];
1848 const IntPolyh_Point &P2=TPoints1[Tri1.SecondPoint()];
1849 const IntPolyh_Point &P3=TPoints1[Tri1.ThirdPoint()];
1850 const IntPolyh_Point &Q1=TPoints2[Tri2.FirstPoint()];
1851 const IntPolyh_Point &Q2=TPoints2[Tri2.SecondPoint()];
1852 const IntPolyh_Point &Q3=TPoints2[Tri2.ThirdPoint()];
1854 /* The first triangle is (p1,p2,p3). The other is (q1,q2,q3).
1855 The edges are (e1,e2,e3) and (f1,f2,f3).
1856 The normals are n1 and m1*/
1858 const IntPolyh_Point e1=P2-P1;
1859 const IntPolyh_Point e2=P3-P2;
1860 const IntPolyh_Point e3=P1-P3;
1862 const IntPolyh_Point f1=Q2-Q1;
1863 const IntPolyh_Point f2=Q3-Q2;
1864 const IntPolyh_Point f3=Q1-Q3;
1866 IntPolyh_Point nn1,mm1;
1867 nn1.Cross(e1, e2); //normal to the first triangle
1868 mm1.Cross(f1, f2); //normal to the second triangle
1870 Standard_Real nn1modulus, mm1modulus;
1871 nn1modulus=sqrt(nn1.SquareModulus());
1872 mm1modulus=sqrt(mm1.SquareModulus());
1874 //-------------------------------------------------
1875 ///calculation of intersections points between triangles
1876 //-------------------------------------------------
1878 Standard_Integer NbPoints=0;
1879 IntPolyh_StartPoint SP1,SP2;
1882 if(Abs(nn1modulus)<MyConfusionPrecision) {//10.0e-20){
1886 const IntPolyh_Point n1=nn1.Divide(nn1modulus);
1888 if( (NbPointsTotal<3)&&(EdgeInit2!=Tri2.FirstEdge()) ) {
1889 IntPolyh_StartPoint Pt1,Pt2;
1890 NbPoints=TriangleEdgeContact2(1,1,Tri1,Tri2,P1,P2,P3,e1,e2,e3,Q1,Q2,f1,n1,Pt1,Pt2);
1891 TestNbPoints(1,NbPoints,NbPointsTotal,Pt1,Pt2,SP1,SP2);
1894 if( (NbPointsTotal<3)&&(EdgeInit2!=Tri2.SecondEdge()) ) {
1895 IntPolyh_StartPoint Pt1,Pt2;
1896 NbPoints=TriangleEdgeContact2(1,2,Tri1,Tri2,P1,P2,P3,e1,e2,e3,Q2,Q3,f2,n1,Pt1,Pt2);
1897 TestNbPoints(1,NbPoints,NbPointsTotal,Pt1,Pt2,SP1,SP2);
1900 if( (NbPointsTotal<3)&&(EdgeInit2!=Tri2.ThirdEdge()) ) {
1901 IntPolyh_StartPoint Pt1,Pt2;
1902 NbPoints=TriangleEdgeContact2(1,3,Tri1,Tri2,P1,P2,P3,e1,e2,e3,Q3,Q1,f3,n1,Pt1,Pt2);
1903 TestNbPoints(1,NbPoints,NbPointsTotal,Pt1,Pt2,SP1,SP2);
1907 if(Abs(mm1modulus)<MyConfusionPrecision) {//10.0e-20){
1911 const IntPolyh_Point m1=mm1.Divide(mm1modulus);
1913 if( (NbPointsTotal<3)&&(EdgeInit1!=Tri1.FirstEdge()) ) {
1914 IntPolyh_StartPoint Pt1,Pt2;
1915 NbPoints=TriangleEdgeContact2(2,1,Tri1,Tri2,Q1,Q2,Q3,f1,f2,f3,P1,P2,e1,m1,Pt1,Pt2);
1916 TestNbPoints(2,NbPoints,NbPointsTotal,Pt1,Pt2,SP1,SP2);
1919 if( (NbPointsTotal<3)&&(EdgeInit1!=Tri1.SecondEdge()) ) {
1920 IntPolyh_StartPoint Pt1,Pt2;
1921 NbPoints=TriangleEdgeContact2(2,2,Tri1,Tri2,Q1,Q2,Q3,f1,f2,f3,P2,P3,e2,m1,Pt1,Pt2);
1922 TestNbPoints(2,NbPoints,NbPointsTotal,Pt1,Pt2,SP1,SP2);
1925 if( (NbPointsTotal<3)&&(EdgeInit1!=Tri1.ThirdEdge()) ) {
1926 IntPolyh_StartPoint Pt1,Pt2;
1927 NbPoints=TriangleEdgeContact2(2,3,Tri1,Tri2,Q1,Q2,Q3,f1,f2,f3,P3,P1,e3,m1,Pt1,Pt2);
1928 TestNbPoints(2,NbPoints,NbPointsTotal,Pt1,Pt2,SP1,SP2);
1932 if (NbPointsTotal==1) {
1933 if(SP1.CheckSameSP(SPInit))
1939 else if( (NbPointsTotal==2)&&(SP1.CheckSameSP(SPInit)) ) {
1940 NbPointsTotal=1;//SP1 et SPInit sont identiques
1943 else if( (NbPointsTotal==2)&&(SP2.CheckSameSP(SPInit)) ) {
1944 NbPointsTotal=1;//SP2 et SPInit sont identiques
1948 else if(NbPointsTotal>1) {
1952 SPNext.SetCoupleValue(T1,T2);
1953 return (NbPointsTotal);
1955 //=======================================================================
1956 //function : CalculPtsInterTriEdgeCoplanaires
1958 //=======================================================================
1959 void CalculPtsInterTriEdgeCoplanaires(const Standard_Integer TriSurfID,
1960 const IntPolyh_Point &NormaleTri,
1961 const IntPolyh_Point &PE1,
1962 const IntPolyh_Point &PE2,
1963 const IntPolyh_Point &Edge,
1964 const IntPolyh_Point &PT1,
1965 const IntPolyh_Point &PT2,
1966 const IntPolyh_Point &Cote,
1967 const Standard_Integer CoteIndex,
1968 IntPolyh_StartPoint &SP1,
1969 IntPolyh_StartPoint &SP2,
1970 Standard_Integer &NbPoints)
1972 IntPolyh_Point TestParalleles;
1973 TestParalleles.Cross(Edge,Cote);
1974 if(sqrt(TestParalleles.SquareModulus())<MyConfusionPrecision) {
1976 Per.Cross(NormaleTri,Cote);
1977 Standard_Real p1p = Per.Dot(PE1);
1978 Standard_Real p2p = Per.Dot(PE2);
1979 Standard_Real p0p = Per.Dot(PT1);
1980 if ( ( (p1p>=p0p)&&(p2p<=p0p) )||( (p1p<=p0p)&&(p2p>=p0p) ) ) {
1981 Standard_Real lambda=(p1p-p0p)/(p1p-p2p);
1982 if (lambda<-MyConfusionPrecision) {
1985 IntPolyh_Point PIE=PE1+Edge*lambda;
1986 Standard_Real alpha=RealLast();
1987 if(Cote.X()!=0) alpha=(PIE.X()-PT1.X())/Cote.X();
1988 else if (Cote.Y()!=0) alpha=(PIE.Y()-PT1.Y())/Cote.Y();
1989 else if (Cote.Z()!=0) alpha=(PIE.Z()-PT1.Z())/Cote.Z();
1993 if (alpha<-MyConfusionPrecision) {
1998 SP1.SetXYZ(PIE.X(),PIE.Y(),PIE.Z());
2000 SP1.SetUV1(PT1.U()+Cote.U()*alpha,PT1.V()+Cote.V()*alpha);
2001 SP1.SetUV2(PIE.U(),PIE.V());
2002 SP1.SetEdge1(CoteIndex);
2005 else if (TriSurfID==2) {
2006 SP1.SetUV1(PIE.U(),PIE.V());
2007 SP1.SetUV2(PT1.U()+Cote.U()*alpha,PT1.V()+Cote.V()*alpha);
2008 SP1.SetEdge2(CoteIndex);
2016 else if (NbPoints==1) {
2017 SP2.SetXYZ(PIE.X(),PIE.Y(),PIE.Z());
2019 SP2.SetUV1(PT1.U()+Cote.U()*alpha,PT1.V()+Cote.V()*alpha);
2020 SP2.SetUV2(PIE.U(),PIE.V());
2021 SP2.SetEdge1(CoteIndex);
2024 else if (TriSurfID==2) {
2025 SP2.SetUV1(PIE.U(),PIE.V());
2026 SP2.SetUV2(PT1.U()+Cote.U()*alpha,PT1.V()+Cote.V()*alpha);
2027 SP2.SetEdge2(CoteIndex);
2035 else if( (NbPoints>2)||(NbPoints<0) ) {
2041 else { //Cote et Edge paralleles, avec les rejections precedentes ils sont sur la meme droite
2042 //On projette les points sur cette droite
2043 Standard_Real pe1p= Cote.Dot(PE1);
2044 Standard_Real pe2p= Cote.Dot(PE2);
2045 Standard_Real pt1p= Cote.Dot(PT1);
2046 Standard_Real pt2p= Cote.Dot(PT2);
2048 IntPolyh_Point PEP1,PTP1,PEP2,PTP2;
2050 //PEP1 et PEP2 sont les points de contact entre le triangle et l'edge dans le repere UV de l'edge
2051 //PTP1 et PTP2 sont les correspondants respectifs a PEP1 et PEP2 dans le repere UV du triangle
2055 if ( (pt1p<pe1p)&&(pe1p<=pt2p) ) {
2057 PTP1=PT1+Cote*((pe1p-pt1p)/(pt2p-pt1p));
2061 PTP2=PT1+Cote*((pe2p-pt1p)/(pt2p-pt1p));
2065 PEP2=PE1+Edge*((pt1p-pe1p)/(pe2p-pe1p));
2070 else if( (pt2p<pe1p)&&(pe1p<=pt1p) ) {
2072 PTP1=PT1+Cote*((pt1p-pe1p)/(pt1p-pt2p));
2076 PTP2=PT1+Cote*((pe2p-pt1p)/(pt2p-pt1p));
2080 PEP2=PE1+Edge*((pt2p-pe1p)/(pe2p-pe1p));
2088 if ( (pt1p<pe2p)&&(pe2p<=pt2p) ) {
2090 PTP1=PT1+Cote*((pe2p-pt1p)/(pt2p-pt1p));
2094 PTP2=PT1+Cote*((pe1p-pt1p)/(pt2p-pt1p));
2098 PEP2=PE2+Edge*((pt1p-pe1p)/(pe2p-pe1p));
2103 else if( (pt2p<pe2p)&&(pe2p<=pt1p) ) {
2105 PTP1=PT1+Cote*((pt1p-pe2p)/(pt1p-pt2p));
2109 PTP2=PT1+Cote*((pe1p-pt1p)/(pt2p-pt1p));
2113 PEP2=PE1+Edge*((pt2p-pe1p)/(pe2p-pe1p));
2121 if (Abs(PEP1.U()-PEP2.U())<MyConfusionPrecision
2122 &&(Abs(PEP1.V()-PEP2.V())<MyConfusionPrecision) ) NbPoints=1;
2124 SP1.SetXYZ(PEP1.X(),PEP1.Y(),PEP1.Z());
2126 SP1.SetUV1(PTP1.U(),PTP1.V());
2127 SP1.SetUV2(PEP1.U(),PEP1.V());
2128 SP1.SetEdge1(CoteIndex);
2131 SP1.SetUV1(PEP1.U(),PTP1.V());
2132 SP1.SetUV2(PTP1.U(),PEP1.V());
2133 SP1.SetEdge2(CoteIndex);
2137 SP2.SetXYZ(PEP2.X(),PEP2.Y(),PEP2.Z());
2139 SP2.SetUV1(PTP2.U(),PTP2.V());
2140 SP2.SetUV2(PEP2.U(),PEP2.V());
2141 SP2.SetEdge1(CoteIndex);
2144 SP2.SetUV1(PEP2.U(),PTP2.V());
2145 SP2.SetUV2(PTP2.U(),PEP2.V());
2146 SP2.SetEdge2(CoteIndex);
2152 //=======================================================================
2153 //function : CalculPtsInterTriEdgeCoplanaires2
2155 //=======================================================================
2156 void CalculPtsInterTriEdgeCoplanaires2(const Standard_Integer TriSurfID,
2157 const IntPolyh_Point &NormaleTri,
2158 const IntPolyh_Triangle &Tri1,
2159 const IntPolyh_Triangle &Tri2,
2160 const IntPolyh_Point &PE1,
2161 const IntPolyh_Point &PE2,
2162 const IntPolyh_Point &Edge,
2163 const Standard_Integer EdgeIndex,
2164 const IntPolyh_Point &PT1,
2165 const IntPolyh_Point &PT2,
2166 const IntPolyh_Point &Cote,
2167 const Standard_Integer CoteIndex,
2168 IntPolyh_StartPoint &SP1,
2169 IntPolyh_StartPoint &SP2,
2170 Standard_Integer &NbPoints)
2172 Standard_Real aDE, aDC;
2174 gp_Vec aVE(Edge.X(), Edge.Y(), Edge.Z());
2175 gp_Vec aVC(Cote.X(), Cote.Y(), Cote.Z());
2177 aDE = aVE.SquareMagnitude();
2178 aDC = aVC.SquareMagnitude();
2180 if (aDE > SquareMyConfusionPrecision) {
2184 if (aDC > SquareMyConfusionPrecision) {
2188 if (!aVE.IsParallel(aVC, MyConfusionPrecision)) {
2189 ///Edge and side are not parallel
2191 Per.Cross(NormaleTri,Cote);
2192 Standard_Real p1p = Per.Dot(PE1);
2193 Standard_Real p2p = Per.Dot(PE2);
2194 Standard_Real p0p = Per.Dot(PT1);
2195 ///The edge are PT1 are projected on the perpendicular of the side in the plane of the triangle
2196 if ( ( (p1p>=p0p)&&(p0p>=p2p) )||( (p1p<=p0p)&&(p0p<=p2p) ) ) {
2197 Standard_Real lambda=(p1p-p0p)/(p1p-p2p);
2198 if (lambda<-MyConfusionPrecision) {
2202 if (Abs(lambda)<MyConfusionPrecision)//lambda=0
2204 else if (Abs(lambda)>1.0-MyConfusionPrecision)//lambda=1
2207 PIE=PE1+Edge*lambda;
2209 Standard_Real alpha=RealLast();
2210 if(Cote.X()!=0) alpha=(PIE.X()-PT1.X())/Cote.X();
2211 else if (Cote.Y()!=0) alpha=(PIE.Y()-PT1.Y())/Cote.Y();
2212 else if (Cote.Z()!=0) alpha=(PIE.Z()-PT1.Z())/Cote.Z();
2217 if (alpha<-MyConfusionPrecision) {
2222 SP1.SetXYZ(PIE.X(),PIE.Y(),PIE.Z());
2224 if(Abs(alpha)<MyConfusionPrecision) {//alpha=0
2225 SP1.SetUV1(PT1.U(),PT1.V());
2226 SP1.SetUV1(PIE.U(),PIE.V());
2229 if(Abs(alpha)>1.0-MyConfusionPrecision) {//alpha=1
2230 SP1.SetUV1(PT2.U(),PT2.V());
2231 SP1.SetUV1(PIE.U(),PIE.V());
2235 SP1.SetUV1(PT1.U()+Cote.U()*alpha,PT1.V()+Cote.V()*alpha);
2236 SP1.SetUV2(PIE.U(),PIE.V());
2237 SP1.SetEdge1(Tri1.GetEdgeNumber(CoteIndex));
2238 if (Tri1.GetEdgeOrientation(CoteIndex)>0) SP1.SetLambda1(alpha);
2239 else SP1.SetLambda1(1.0-alpha);
2243 else if (TriSurfID==2) {
2244 if(Abs(alpha)<MyConfusionPrecision) {//alpha=0
2245 SP1.SetUV1(PT1.U(),PT1.V());
2246 SP1.SetUV1(PIE.U(),PIE.V());
2249 if(Abs(alpha)>1.0-MyConfusionPrecision) {//alpha=1
2250 SP1.SetUV1(PT2.U(),PT2.V());
2251 SP1.SetUV1(PIE.U(),PIE.V());
2255 SP1.SetUV1(PIE.U(),PIE.V());
2256 SP1.SetUV2(PT1.U()+Cote.U()*alpha,PT1.V()+Cote.V()*alpha);
2257 SP1.SetEdge2(Tri2.GetEdgeNumber(CoteIndex));
2258 if (Tri2.GetEdgeOrientation(CoteIndex)>0) SP1.SetLambda2(alpha);
2259 else SP1.SetLambda2(1.0-alpha);
2268 else if (NbPoints==1) {
2269 SP2.SetXYZ(PIE.X(),PIE.Y(),PIE.Z());
2271 if(Abs(alpha)<MyConfusionPrecision) {//alpha=0
2272 SP2.SetUV1(PT1.U(),PT1.V());
2273 SP2.SetUV1(PIE.U(),PIE.V());
2276 if(Abs(alpha)>1.0-MyConfusionPrecision) {//alpha=1
2277 SP2.SetUV1(PT2.U(),PT2.V());
2278 SP2.SetUV1(PIE.U(),PIE.V());
2282 SP2.SetUV1(PT1.U()+Cote.U()*alpha,PT1.V()+Cote.V()*alpha);
2283 SP2.SetUV2(PIE.U(),PIE.V());
2284 SP2.SetEdge1(Tri1.GetEdgeNumber(CoteIndex));
2285 if (Tri1.GetEdgeOrientation(CoteIndex)>0) SP2.SetLambda1(alpha);
2286 else SP2.SetLambda1(1.0-alpha);
2290 else if (TriSurfID==2) {
2291 if(Abs(alpha)<MyConfusionPrecision) {//alpha=0
2292 SP2.SetUV1(PT1.U(),PT1.V());
2293 SP2.SetUV1(PIE.U(),PIE.V());
2296 if(Abs(alpha)>1.0-MyConfusionPrecision) {//alpha=1
2297 SP2.SetUV1(PT2.U(),PT2.V());
2298 SP2.SetUV1(PIE.U(),PIE.V());
2302 SP2.SetUV1(PIE.U(),PIE.V());
2303 SP2.SetUV2(PT1.U()+Cote.U()*alpha,PT1.V()+Cote.V()*alpha);
2304 SP2.SetEdge2(Tri2.GetEdgeNumber(CoteIndex));
2305 if (Tri1.GetEdgeOrientation(CoteIndex)>0) SP2.SetLambda2(alpha);
2306 else SP2.SetLambda2(1.0-alpha);
2315 else if( (NbPoints>2)||(NbPoints<0) ) {
2322 //Side and Edge are parallel, with previous
2323 //rejections they are at the same side
2324 //The points are projected on that side
2325 Standard_Real pe1p= Cote.Dot(PE1);
2326 Standard_Real pe2p= Cote.Dot(PE2);
2327 Standard_Real pt1p= Cote.Dot(PT1);
2328 Standard_Real pt2p= Cote.Dot(PT2);
2329 Standard_Real lambda1=0., lambda2=0., alpha1=0., alpha2=0.;
2330 IntPolyh_Point PEP1,PTP1,PEP2,PTP2;
2333 if ( (pt1p<pe1p)&&(pe1p<=pt2p) ) {
2336 alpha1=((pe1p-pt1p)/(pt2p-pt1p));
2337 PTP1=PT1+Cote*alpha1;
2342 alpha2=((pe2p-pt1p)/(pt2p-pt1p));
2343 PTP2=PT1+Cote*alpha2;
2347 lambda2=((pt1p-pe1p)/(pe2p-pe1p));
2348 PEP2=PE1+Edge*lambda2;
2354 else if( (pt2p<pe1p)&&(pe1p<=pt1p) ) {
2357 alpha1=((pt1p-pe1p)/(pt1p-pt2p));
2358 PTP1=PT1+Cote*alpha1;
2363 alpha2=((pe2p-pt1p)/(pt2p-pt1p));
2364 PTP2=PT1+Cote*alpha2;
2368 lambda2=((pt2p-pe1p)/(pe2p-pe1p));
2369 PEP2=PE1+Edge*lambda2;
2378 if ( (pt1p<pe2p)&&(pe2p<=pt2p) ) {
2381 alpha1=((pe2p-pt1p)/(pt2p-pt1p));
2382 PTP1=PT1+Cote*alpha1;
2387 alpha2=((pe1p-pt1p)/(pt2p-pt1p));
2388 PTP2=PT1+Cote*alpha2;
2392 lambda2=((pt1p-pe1p)/(pe2p-pe1p));
2393 PEP2=PE2+Edge*lambda2;
2399 else if( (pt2p<pe2p)&&(pe2p<=pt1p) ) {
2402 alpha1=((pt1p-pe2p)/(pt1p-pt2p));
2403 PTP1=PT1+Cote*alpha1;
2408 alpha2=((pe1p-pt1p)/(pt2p-pt1p));
2409 PTP2=PT1+Cote*alpha2;
2413 lambda2=((pt2p-pe1p)/(pe2p-pe1p));
2414 PEP2=PE1+Edge*lambda2;
2423 SP1.SetXYZ(PEP1.X(),PEP1.Y(),PEP1.Z());
2424 if (TriSurfID==1) {///cote appartient a Tri1
2425 SP1.SetUV1(PTP1.U(),PTP1.V());
2426 SP1.SetUV2(PEP1.U(),PEP1.V());
2427 SP1.SetEdge1(Tri1.GetEdgeNumber(CoteIndex));
2429 if(Tri1.GetEdgeOrientation(CoteIndex)>0) SP1.SetLambda1(alpha1);
2430 else SP1.SetLambda1(1.0-alpha1);
2432 if(Tri2.GetEdgeOrientation(EdgeIndex)>0) SP1.SetLambda2(lambda1);
2433 else SP1.SetLambda2(1.0-lambda1);
2435 if (TriSurfID==2) {///cote appartient a Tri2
2436 SP1.SetUV1(PEP1.U(),PTP1.V());
2437 SP1.SetUV2(PTP1.U(),PEP1.V());
2438 SP1.SetEdge2(Tri2.GetEdgeNumber(CoteIndex));
2440 if(Tri2.GetEdgeOrientation(CoteIndex)>0) SP1.SetLambda1(alpha1);
2441 else SP1.SetLambda1(1.0-alpha1);
2443 if(Tri1.GetEdgeOrientation(EdgeIndex)>0) SP1.SetLambda2(lambda1);
2444 else SP1.SetLambda2(1.0-lambda1);
2447 //It is checked if PEP1!=PEP2
2448 if ( (NbPoints==2)&&(Abs(PEP1.U()-PEP2.U())<MyConfusionPrecision)
2449 &&(Abs(PEP1.V()-PEP2.V())<MyConfusionPrecision) ) NbPoints=1;
2451 SP2.SetXYZ(PEP2.X(),PEP2.Y(),PEP2.Z());
2453 SP2.SetUV1(PTP2.U(),PTP2.V());
2454 SP2.SetUV2(PEP2.U(),PEP2.V());
2455 SP2.SetEdge1(Tri1.GetEdgeNumber(CoteIndex));
2457 if(Tri1.GetEdgeOrientation(CoteIndex)>0) SP2.SetLambda1(alpha1);
2458 else SP2.SetLambda1(1.0-alpha1);
2460 if(Tri2.GetEdgeOrientation(EdgeIndex)>0) SP2.SetLambda2(lambda1);
2461 else SP2.SetLambda2(1.0-lambda1);
2464 SP2.SetUV1(PEP2.U(),PTP2.V());
2465 SP2.SetUV2(PTP2.U(),PEP2.V());
2466 SP2.SetEdge2(Tri2.GetEdgeNumber(CoteIndex));
2468 if(Tri1.GetEdgeOrientation(CoteIndex)>0) SP2.SetLambda1(alpha1);
2469 else SP2.SetLambda1(1.0-alpha1);
2471 if(Tri2.GetEdgeOrientation(EdgeIndex)>0) SP2.SetLambda2(lambda1);
2472 else SP2.SetLambda2(1.0-lambda1);
2477 //Filter if the point is placed on top, the edge is set to -1
2479 if(Abs(SP1.Lambda1())<MyConfusionPrecision)
2481 if(Abs(SP1.Lambda1()-1)<MyConfusionPrecision)
2483 if(Abs(SP1.Lambda2())<MyConfusionPrecision)
2485 if(Abs(SP1.Lambda2()-1)<MyConfusionPrecision)
2489 if(Abs(SP2.Lambda1())<MyConfusionPrecision)
2491 if(Abs(SP2.Lambda1()-1)<MyConfusionPrecision)
2493 if(Abs(SP2.Lambda2())<MyConfusionPrecision)
2495 if(Abs(SP2.Lambda2()-1)<MyConfusionPrecision)
2499 //=======================================================================
2500 //function : TriangleEdgeContact
2502 //=======================================================================
2503 Standard_Integer IntPolyh_MaillageAffinage::TriangleEdgeContact
2504 (const Standard_Integer TriSurfID,
2505 const Standard_Integer EdgeIndex,
2506 const IntPolyh_Point &PT1,
2507 const IntPolyh_Point &PT2,
2508 const IntPolyh_Point &PT3,
2509 const IntPolyh_Point &Cote12,
2510 const IntPolyh_Point &Cote23,
2511 const IntPolyh_Point &Cote31,
2512 const IntPolyh_Point &PE1,const IntPolyh_Point &PE2,
2513 const IntPolyh_Point &Edge,
2514 const IntPolyh_Point &NormaleT,
2515 IntPolyh_StartPoint &SP1,IntPolyh_StartPoint &SP2) const
2518 Standard_Real lambda =0.;
2519 Standard_Real alpha =0.;
2520 Standard_Real beta =0.;
2523 //The edge, on which the point is located, is known.
2525 SP1.SetEdge2(EdgeIndex);
2526 SP2.SetEdge2(EdgeIndex);
2528 else if (TriSurfID==2) {
2529 SP1.SetEdge1(EdgeIndex);
2530 SP2.SetEdge1(EdgeIndex);
2536 /**The edge is projected on the normal of the triangle if yes
2537 --> free intersection (point I)--> start point is found*/
2538 Standard_Integer NbPoints=0;
2540 if(NormaleT.SquareModulus()==0) {
2543 else if( (Cote12.SquareModulus()==0)
2544 ||(Cote23.SquareModulus()==0)
2545 ||(Cote31.SquareModulus()==0) ) {
2548 else if(Edge.SquareModulus()==0) {
2552 Standard_Real pe1 = NormaleT.Dot(PE1);
2553 Standard_Real pe2 = NormaleT.Dot(PE2);
2554 Standard_Real pt1 = NormaleT.Dot(PT1);
2556 // PE1I = lambda.Edge
2558 if( (Abs(pe1-pe2)<MyConfusionPrecision)&&(Abs(pe1-pt1)<MyConfusionPrecision) ) {
2559 //edge and triangle are coplanar (two contact points maximum)
2561 //The tops of the triangle are projected on the perpendicular of the edge
2563 IntPolyh_Point PerpEdge;
2564 PerpEdge.Cross(NormaleT,Edge);
2565 Standard_Real pp1 = PerpEdge.Dot(PT1);
2566 Standard_Real pp2 = PerpEdge.Dot(PT2);
2567 Standard_Real pp3 = PerpEdge.Dot(PT3);
2568 Standard_Real ppe1 = PerpEdge.Dot(PE1);
2570 if ( ( (pp1>ppe1)&&(pp2<=ppe1)&&(pp3<=ppe1) ) || ( (pp1<ppe1)&&(pp2>=ppe1)&&(pp3>=ppe1) ) ){
2571 //there are two sides (commun top PT1) that can cut the edge
2574 CalculPtsInterTriEdgeCoplanaires(TriSurfID,NormaleT,
2575 PE1,PE2,Edge,PT1,PT2,
2576 Cote12,1,SP1,SP2,NbPoints);
2578 if ( (NbPoints>1)&&(Abs(SP2.U1()-SP1.U1())<MyConfusionPrecision)
2579 &&(Abs(SP1.V1()-SP2.V1())<MyConfusionPrecision) ) NbPoints=1;
2583 CalculPtsInterTriEdgeCoplanaires(TriSurfID,NormaleT,
2584 PE1,PE2,Edge,PT3,PT1,
2585 Cote31,3,SP1,SP2,NbPoints);
2589 if ( (NbPoints>1)&&(Abs(SP1.U1()-SP2.U1())<MyConfusionPrecision)
2590 &&(Abs(SP1.V2()-SP2.V1())<MyConfusionPrecision) ) NbPoints=1;
2591 if (NbPoints>=2) return(NbPoints);
2593 else if ( ( ( (pp2>ppe1)&&(pp1<=ppe1)&&(pp3<=ppe1) ) || ( (pp2<ppe1)&&(pp1>=ppe1)&&(pp3>=ppe1) ) )
2595 //there are two sides (common top PT2) that can cut the edge
2598 CalculPtsInterTriEdgeCoplanaires(TriSurfID,NormaleT,
2599 PE1,PE2,Edge,PT1,PT2,
2600 Cote12,1,SP1,SP2,NbPoints);
2602 if ( (NbPoints>1)&&(Abs(SP2.U1()-SP1.U1())<MyConfusionPrecision)
2603 &&(Abs(SP1.V1()-SP2.V1())<MyConfusionPrecision) ) NbPoints=1;
2606 if(NbPoints<2) CalculPtsInterTriEdgeCoplanaires(TriSurfID,NormaleT,
2607 PE1,PE2,Edge,PT2,PT3,
2608 Cote23,2,SP1,SP2,NbPoints);
2610 if ( (NbPoints>1)&&(Abs(SP2.U1()-SP1.U1())<MyConfusionPrecision)
2611 &&(Abs(SP1.V1()-SP2.V1())<MyConfusionPrecision) ) NbPoints=1;
2612 if (NbPoints>=2) return(NbPoints);
2614 else if ( (( (pp3>ppe1)&&(pp1<=ppe1)&&(pp2<=ppe1) ) || ( (pp3<ppe1)&&(pp1>=ppe1)&&(pp2>=ppe1) ))
2616 //there are two sides (common top PT3) that can cut the edge
2619 CalculPtsInterTriEdgeCoplanaires(TriSurfID,NormaleT,
2620 PE1,PE2,Edge,PT3,PT1,Cote31,
2621 3,SP1,SP2,NbPoints);
2623 if ( (NbPoints>1)&&(Abs(SP2.U1()-SP1.U1())<MyConfusionPrecision)
2624 &&(Abs(SP1.V1()-SP2.V1())<MyConfusionPrecision) ) NbPoints=1;
2627 if (NbPoints<2) CalculPtsInterTriEdgeCoplanaires(TriSurfID,NormaleT,
2628 PE1,PE2,Edge,PT2,PT3,
2629 Cote23,2,SP1,SP2,NbPoints);
2631 if ( (NbPoints>1)&&(Abs(SP2.U1()-SP1.U1())<MyConfusionPrecision)
2632 &&(Abs(SP2.V1()-SP1.V1())<MyConfusionPrecision) ) NbPoints=1;
2633 if (NbPoints>=2) return(NbPoints);
2637 //------------------------------------------------------
2638 // edge and triangle NON COPLANAR (a contact point)
2639 //------------------------------------------------------
2640 else if( ( (pe1>=pt1)&&(pe2<=pt1) ) || ( (pe1<=pt1)&&(pe2>=pt1) ) ) {
2641 lambda=(pe1-pt1)/(pe1-pe2);
2643 if (lambda<-MyConfusionPrecision) {
2646 else if (Abs(lambda)<MyConfusionPrecision) {//lambda==0
2648 if(TriSurfID==1) SP1.SetEdge2(0);
2649 else SP1.SetEdge1(0);
2651 else if (Abs(lambda-1.0)<MyConfusionPrecision) {//lambda==1
2653 if(TriSurfID==1) SP1.SetEdge2(0);
2654 else SP1.SetEdge1(0);
2658 if(TriSurfID==1) SP1.SetEdge2(EdgeIndex);
2659 else SP1.SetEdge1(EdgeIndex);
2662 if(Abs(Cote23.X())>MyConfusionPrecision) {
2663 Standard_Real D=(Cote12.Y()-Cote12.X()*Cote23.Y()/Cote23.X());
2664 if(D!=0) alpha = (PI.Y()-PT1.Y()-(PI.X()-PT1.X())*Cote23.Y()/Cote23.X())/D;
2668 if ((alpha<-MyConfusionPrecision)||(alpha>(1.0+MyConfusionPrecision))) return(0);
2669 else beta = (PI.X()-PT1.X()-alpha*Cote12.X())/Cote23.X();
2672 else if (Abs(Cote12.X())>MyConfusionPrecision) { //On a Cote23.X()==0
2673 alpha = (PI.X()-PT1.X())/Cote12.X();
2675 if ((alpha<-MyConfusionPrecision)||(alpha>(1.0+MyConfusionPrecision))) return(0);
2677 else if (Abs(Cote23.Y())>MyConfusionPrecision) beta = (PI.Y()-PT1.Y()-alpha*Cote12.Y())/Cote23.Y();
2678 else if (Abs(Cote23.Z())>MyConfusionPrecision) beta = (PI.Z()-PT1.Z()-alpha*Cote12.Z())/Cote23.Z();
2684 else if (Abs(Cote23.Y())>MyConfusionPrecision) {
2685 //On a Cote23.X()==0 et Cote12.X()==0 ==> equation can't be used
2686 Standard_Real D=(Cote12.Z()-Cote12.Y()*Cote23.Z()/Cote23.Y());
2688 if(D!=0) alpha = (PI.Z()-PT1.Z()-(PI.Y()-PT1.Y())*Cote23.Z()/Cote23.Y())/D;
2693 if ( (alpha<-MyConfusionPrecision)||(alpha>(1.0+MyConfusionPrecision)) ) return(0);
2694 else beta = (PI.Y()-PT1.Y()-alpha*Cote12.Y())/Cote23.Y();
2697 else if (Abs(Cote12.Y())>MyConfusionPrecision) {
2698 //On a Cote23.X()==0, Cote12.X()==0 et Cote23.Y()==0
2699 alpha = (PI.Y()-PT1.Y())/Cote12.Y();
2701 if ((Abs(alpha)<MyConfusionPrecision)||(Abs(alpha-1.0)<MyConfusionPrecision)) return(0);
2703 else if (Abs(Cote23.Z())>MyConfusionPrecision) beta = (PI.Z()-PT1.Z()-alpha*Cote12.Z())/Cote23.Z();
2710 else { //two equations of three can't be used
2716 if( (beta<-MyConfusionPrecision)||(beta>(alpha+MyConfusionPrecision)) ) return(0);
2718 SP1.SetXYZ(PI.X(),PI.Y(),PI.Z());
2721 SP1.SetUV2(PI.U(),PI.V());
2722 SP1.SetUV1(PT1.U()+Cote12.U()*alpha+Cote23.U()*beta, PT1.V()+Cote12.V()*alpha+Cote23.V()*beta);
2724 if (beta<MyConfusionPrecision) {//beta==0 && alpha
2726 SP1.SetLambda1(alpha);
2728 if (Abs(beta-alpha)<MyConfusionPrecision) {//beta==alpha
2730 SP1.SetLambda1(1.0-alpha);
2732 if (Abs(alpha-1)<MyConfusionPrecision)//alpha==1
2734 if (alpha<MyConfusionPrecision) {//alpha=0 --> beta==0
2735 SP1.SetXYZ(PT1.X(),PT1.Y(),PT1.Z());
2736 SP1.SetUV1(PT1.U(),PT1.V());
2739 if ( (beta<MyConfusionPrecision)&&(Abs(1-alpha)<MyConfusionPrecision) ) {//beta==0 alpha==1
2740 SP1.SetXYZ(PT2.X(),PT2.Y(),PT2.Z());
2741 SP1.SetUV1(PT2.U(),PT2.V());
2744 if ( (Abs(beta-1)<MyConfusionPrecision)||(Abs(1-alpha)<MyConfusionPrecision) ) {//beta==1 alpha==1
2745 SP1.SetXYZ(PT3.X(),PT3.Y(),PT3.Z());
2746 SP1.SetUV1(PT3.U(),PT3.V());
2750 else if(TriSurfID==2) {
2751 SP1.SetUV1(PI.U(),PI.V());
2752 SP1.SetUV2(PT1.U()+Cote12.U()*alpha+Cote23.U()*beta, PT1.V()+Cote12.V()*alpha+Cote23.V()*beta);
2754 if (beta<MyConfusionPrecision) { //beta==0
2757 if (Abs(beta-alpha)<MyConfusionPrecision)//beta==alpha
2759 if (Abs(alpha-1)<MyConfusionPrecision)//alpha==1
2761 if (alpha<MyConfusionPrecision) {//alpha=0 --> beta==0
2762 SP1.SetXYZ(PT1.X(),PT1.Y(),PT1.Z());
2763 SP1.SetUV2(PT1.U(),PT1.V());
2766 if ( (beta<MyConfusionPrecision)&&(Abs(1-alpha)<MyConfusionPrecision) ) {//beta==0 alpha==1
2767 SP1.SetXYZ(PT2.X(),PT2.Y(),PT2.Z());
2768 SP1.SetUV2(PT2.U(),PT2.V());
2771 if ( (Abs(beta-1)<MyConfusionPrecision)||(Abs(1-alpha)<MyConfusionPrecision) ) {//beta==1 alpha==1
2772 SP1.SetXYZ(PT3.X(),PT3.Y(),PT3.Z());
2773 SP1.SetUV2(PT3.U(),PT3.V());
2787 //=======================================================================
2788 //function : TriangleEdgeContact2
2790 //=======================================================================
2791 Standard_Integer IntPolyh_MaillageAffinage::TriangleEdgeContact2
2792 (const Standard_Integer TriSurfID,
2793 const Standard_Integer EdgeIndex,
2794 const IntPolyh_Triangle &Tri1,
2795 const IntPolyh_Triangle &Tri2,
2796 const IntPolyh_Point &PT1,
2797 const IntPolyh_Point &PT2,
2798 const IntPolyh_Point &PT3,
2799 const IntPolyh_Point &Cote12,
2800 const IntPolyh_Point &Cote23,
2801 const IntPolyh_Point &Cote31,
2802 const IntPolyh_Point &PE1,const IntPolyh_Point &PE2,
2803 const IntPolyh_Point &Edge,
2804 const IntPolyh_Point &NormaleT,
2805 IntPolyh_StartPoint &SP1,IntPolyh_StartPoint &SP2) const
2808 Standard_Real lambda =0., alpha =0., beta =0.;
2810 //One of edges on which the point is located is known
2812 SP1.SetEdge2(Tri2.GetEdgeNumber(EdgeIndex));
2813 SP2.SetEdge2(Tri2.GetEdgeNumber(EdgeIndex));
2815 else if (TriSurfID==2) {
2816 SP1.SetEdge1(Tri1.GetEdgeNumber(EdgeIndex));
2817 SP2.SetEdge1(Tri1.GetEdgeNumber(EdgeIndex));
2823 /**The edge is projected on the normal in the triangle if yes
2824 --> a free intersection (point I)--> Start point is found */
2825 Standard_Integer NbPoints=0;
2826 if(NormaleT.SquareModulus()==0) {
2829 else if( (Cote12.SquareModulus()==0)
2830 ||(Cote23.SquareModulus()==0)
2831 ||(Cote31.SquareModulus()==0) ) {
2834 else if(Edge.SquareModulus()==0) {
2838 Standard_Real pe1 = NormaleT.Dot(PE1);
2839 Standard_Real pe2 = NormaleT.Dot(PE2);
2840 Standard_Real pt1 = NormaleT.Dot(PT1);
2842 // PE1I = lambda.Edge
2844 if( (Abs(pe1-pt1)<MyConfusionPrecision)&&(Abs(pe2-pt1)<MyConfusionPrecision)) {
2845 //edge and triangle are coplanar (two contact points at maximum)
2848 //the tops of the triangle are projected on the perpendicular to the edge
2849 IntPolyh_Point PerpEdge;
2850 PerpEdge.Cross(NormaleT,Edge);
2851 Standard_Real pp1 = PerpEdge.Dot(PT1);
2852 Standard_Real pp2 = PerpEdge.Dot(PT2);
2853 Standard_Real pp3 = PerpEdge.Dot(PT3);
2854 Standard_Real ppe1 = PerpEdge.Dot(PE1);
2857 if ( (Abs(pp1-pp2)<MyConfusionPrecision)&&(Abs(pp1-pp3)<MyConfusionPrecision) ) {
2861 if ( ( (pp1>=ppe1)&&(pp2<=ppe1)&&(pp3<=ppe1) ) || ( (pp1<=ppe1)&&(pp2>=ppe1)&&(pp3>=ppe1) ) ){
2862 //there are two sides (common top PT1) that can cut the edge
2865 CalculPtsInterTriEdgeCoplanaires2(TriSurfID,NormaleT,Tri1,Tri2,PE1,PE2,Edge,EdgeIndex,
2866 PT1,PT2,Cote12,1,SP1,SP2,NbPoints);
2868 if ( (NbPoints>1)&&(Abs(SP2.U1()-SP1.U1())<MyConfusionPrecision)
2869 &&(Abs(SP1.V1()-SP2.V1())<MyConfusionPrecision) ) NbPoints=1;
2872 if (NbPoints<2) CalculPtsInterTriEdgeCoplanaires2(TriSurfID,NormaleT,Tri1,Tri2,PE1,PE2,Edge,EdgeIndex,
2873 PT3,PT1,Cote31,3,SP1,SP2,NbPoints);
2876 if ( (NbPoints>1)&&(Abs(SP1.U1()-SP2.U1())<MyConfusionPrecision)
2877 &&(Abs(SP1.V2()-SP2.V1())<MyConfusionPrecision) ) NbPoints=1;
2878 if (NbPoints>=2) return(NbPoints);
2880 else if ( ( ( (pp2>=ppe1)&&(pp1<=ppe1)&&(pp3<=ppe1) ) || ( (pp2<=ppe1)&&(pp1>=ppe1)&&(pp3>=ppe1) ) )
2882 //there are two sides (common top PT2) that can cut the edge
2885 CalculPtsInterTriEdgeCoplanaires2(TriSurfID,NormaleT,Tri1,Tri2,PE1,PE2,Edge,EdgeIndex,
2886 PT1,PT2,Cote12,1,SP1,SP2,NbPoints);
2888 if ( (NbPoints>1)&&(Abs(SP2.U1()-SP1.U1())<MyConfusionPrecision)
2889 &&(Abs(SP1.V1()-SP2.V1())<MyConfusionPrecision) ) NbPoints=1;
2892 if(NbPoints<2) CalculPtsInterTriEdgeCoplanaires2(TriSurfID,NormaleT,Tri1,Tri2,PE1,PE2,Edge,EdgeIndex,
2893 PT2,PT3,Cote23,2,SP1,SP2,NbPoints);
2895 if ( (NbPoints>1)&&(Abs(SP2.U1()-SP1.U1())<MyConfusionPrecision)
2896 &&(Abs(SP1.V1()-SP2.V1())<MyConfusionPrecision) ) NbPoints=1;
2897 if (NbPoints>=2) return(NbPoints);
2899 else if ( (( (pp3>=ppe1)&&(pp1<=ppe1)&&(pp2<=ppe1) ) || ( (pp3<=ppe1)&&(pp1>=ppe1)&&(pp2>=ppe1) ))
2901 //there are two sides (common top PT3) that can cut the edge
2904 CalculPtsInterTriEdgeCoplanaires2(TriSurfID,NormaleT,Tri1,Tri2,PE1,PE2,Edge,EdgeIndex,
2905 PT3,PT1,Cote31,3,SP1,SP2,NbPoints);
2907 if ( (NbPoints>1)&&(Abs(SP2.U1()-SP1.U1())<MyConfusionPrecision)
2908 &&(Abs(SP1.V1()-SP2.V1())<MyConfusionPrecision) ) NbPoints=1;
2911 if (NbPoints<2) CalculPtsInterTriEdgeCoplanaires2(TriSurfID,NormaleT,Tri1,Tri2,PE1,PE2,Edge,EdgeIndex,
2912 PT2,PT3,Cote23,2,SP1,SP2,NbPoints);
2914 if ( (NbPoints>1)&&(Abs(SP2.U1()-SP1.U1())<MyConfusionPrecision)
2915 &&(Abs(SP2.V1()-SP1.V1())<MyConfusionPrecision) ) NbPoints=1;
2916 if (NbPoints>=2) return(NbPoints);
2920 //------------------------------------------------------
2921 // NON COPLANAR edge and triangle (a contact point)
2922 //------------------------------------------------------
2923 else if( ( (pe1>=pt1)&&(pt1>=pe2) ) || ( (pe1<=pt1)&&(pt1<=pe2) ) ) { //
2924 lambda=(pe1-pt1)/(pe1-pe2);
2926 if (lambda<-MyConfusionPrecision) {
2929 else if (Abs(lambda)<MyConfusionPrecision) {//lambda==0
2931 if(TriSurfID==1) SP1.SetEdge2(-1);
2932 else SP1.SetEdge1(-1);
2934 else if (Abs(lambda-1.0)<MyConfusionPrecision) {//lambda==1
2936 if(TriSurfID==1) SP1.SetEdge2(-1);
2937 else SP1.SetEdge1(-1);
2942 if(Tri2.GetEdgeOrientation(EdgeIndex)>0)
2943 SP1.SetLambda2(lambda);
2944 else SP1.SetLambda2(1.0-lambda);
2947 if(Tri1.GetEdgeOrientation(EdgeIndex)>0)
2948 SP1.SetLambda1(lambda);
2949 else SP1.SetLambda1(1.0-lambda);
2953 Standard_Real Cote23X=Cote23.X();
2954 Standard_Real D1=0.0;
2955 Standard_Real D3,D4;
2957 //Combination Eq1 Eq2
2958 if(Abs(Cote23X)>MyConfusionPrecision) {
2959 D1=Cote12.Y()-Cote12.X()*Cote23.Y()/Cote23X;
2961 if(Abs(D1)>MyConfusionPrecision) {
2962 alpha = ( PI.Y()-PT1.Y()-(PI.X()-PT1.X())*Cote23.Y()/Cote23X )/D1;
2964 ///It is checked if 1.0>=alpha>=0.0
2965 if ((alpha<-MyConfusionPrecision)||(alpha>(1.0+MyConfusionPrecision))) return(0);
2966 else beta = (PI.X()-PT1.X()-alpha*Cote12.X())/Cote23X;
2968 //Combination Eq1 and Eq2 with Cote23.X()==0
2969 else if ( (Abs(Cote12.X())>MyConfusionPrecision)
2970 &&(Abs(Cote23X)<MyConfusionPrecision) ) { //There is Cote23.X()==0
2971 alpha = (PI.X()-PT1.X())/Cote12.X();
2973 if ((alpha<-MyConfusionPrecision)||(alpha>(1.0+MyConfusionPrecision))) return(0);
2975 else if (Abs(Cote23.Y())>MyConfusionPrecision) beta = (PI.Y()-PT1.Y()-alpha*Cote12.Y())/Cote23.Y();
2976 else if (Abs(Cote23.Z())>MyConfusionPrecision) beta = (PI.Z()-PT1.Z()-alpha*Cote12.Z())/Cote23.Z();
2981 //Combination Eq1 and Eq3
2982 else if ( (Abs(Cote23.X())>MyConfusionPrecision)
2983 &&(Abs( D3= (Cote12.Z()-Cote12.X()*Cote23.Z()/Cote23.X()) ) > MyConfusionPrecision) ) {
2985 alpha = (PI.Z()-PT1.Z()-(PI.X()-PT1.X())*Cote23.Z()/Cote23.X())/D3;
2987 if ( (alpha<-MyConfusionPrecision)||(alpha>(1.0+MyConfusionPrecision)) ) return(0);
2988 else beta = (PI.X()-PT1.X()-alpha*Cote12.X())/Cote23.X();
2990 //Combination Eq2 and Eq3
2991 else if ( (Abs(Cote23.Y())>MyConfusionPrecision)
2992 &&(Abs( D4= (Cote12.Z()-Cote12.Y()*Cote23.Z()/Cote23.Y()) ) > MyConfusionPrecision) ) {
2994 alpha = (PI.Z()-PT1.Z()-(PI.Y()-PT1.Y())*Cote23.Z()/Cote23.Y())/D4;
2996 if ( (alpha<-MyConfusionPrecision)||(alpha>(1.0+MyConfusionPrecision)) ) return(0);
2997 else beta = (PI.Y()-PT1.Y()-alpha*Cote12.Y())/Cote23.Y();
2999 //Combination Eq2 and Eq3 with Cote23.Y()==0
3000 else if ( (Abs(Cote12.Y())>MyConfusionPrecision)
3001 && (Abs(Cote23.Y())<MyConfusionPrecision) ) {
3002 alpha = (PI.Y()-PT1.Y())/Cote12.Y();
3004 if ( (alpha<-MyConfusionPrecision)||(alpha>(1.0+MyConfusionPrecision)) ) return(0);
3006 else if (Abs(Cote23.Z())>MyConfusionPrecision) beta = (PI.Z()-PT1.Z()-alpha*Cote12.Z())/Cote23.Z();
3009 printf("\nCote PT2PT3 nul1\n");
3014 //Combination Eq1 and Eq3 with Cote23.Z()==0
3015 else if ( (Abs(Cote12.Z())>MyConfusionPrecision)
3016 && (Abs(Cote23.Z())<MyConfusionPrecision) ) {
3017 alpha = (PI.Z()-PT1.Z())/Cote12.Z();
3019 if ( (alpha<-MyConfusionPrecision)||(alpha>(1.0+MyConfusionPrecision)) ) return(0);
3021 else if (Abs(Cote23.X())>MyConfusionPrecision) beta = (PI.X()-PT1.X()-alpha*Cote12.X())/Cote23.X();
3028 else { //Particular case not processed ?
3034 if( (beta<-MyConfusionPrecision)||(beta>(alpha+MyConfusionPrecision)) ) return(0);
3036 SP1.SetXYZ(PI.X(),PI.Y(),PI.Z());
3039 SP1.SetUV2(PI.U(),PI.V());
3040 SP1.SetUV1(PT1.U()+Cote12.U()*alpha+Cote23.U()*beta, PT1.V()+Cote12.V()*alpha+Cote23.V()*beta);
3042 if (alpha<MyConfusionPrecision) {//alpha=0 --> beta==0
3043 SP1.SetXYZ(PT1.X(),PT1.Y(),PT1.Z());
3044 SP1.SetUV1(PT1.U(),PT1.V());
3047 else if ( (beta<MyConfusionPrecision)&&(Abs(1-alpha)<MyConfusionPrecision) ) {//beta==0 alpha==1
3048 SP1.SetXYZ(PT2.X(),PT2.Y(),PT2.Z());
3049 SP1.SetUV1(PT2.U(),PT2.V());
3052 else if ( (Abs(beta-1)<MyConfusionPrecision)&&(Abs(1-alpha)<MyConfusionPrecision) ) {//beta==1 alpha==1
3053 SP1.SetXYZ(PT3.X(),PT3.Y(),PT3.Z());
3054 SP1.SetUV1(PT3.U(),PT3.V());
3057 else if (beta<MyConfusionPrecision) {//beta==0
3058 SP1.SetEdge1(Tri1.GetEdgeNumber(1));
3059 if(Tri1.GetEdgeOrientation(1)>0)
3060 SP1.SetLambda1(alpha);
3061 else SP1.SetLambda1(1.0-alpha);
3063 else if (Abs(beta-alpha)<MyConfusionPrecision) {//beta==alpha
3064 SP1.SetEdge1(Tri1.GetEdgeNumber(3));
3065 if(Tri1.GetEdgeOrientation(3)>0)
3066 SP1.SetLambda1(1.0-alpha);
3067 else SP1.SetLambda1(alpha);
3069 else if (Abs(alpha-1)<MyConfusionPrecision) {//alpha==1
3070 SP1.SetEdge1(Tri1.GetEdgeNumber(2));
3071 if(Tri1.GetEdgeOrientation(2)>0)
3072 SP1.SetLambda1(beta);
3073 else SP1.SetLambda1(1.0-beta);
3076 else if(TriSurfID==2) {
3077 SP1.SetUV1(PI.U(),PI.V());
3078 SP1.SetUV2(PT1.U()+Cote12.U()*alpha+Cote23.U()*beta, PT1.V()+Cote12.V()*alpha+Cote23.V()*beta);
3080 if (alpha<MyConfusionPrecision) {//alpha=0 --> beta==0
3081 SP1.SetXYZ(PT1.X(),PT1.Y(),PT1.Z());
3082 SP1.SetUV2(PT1.U(),PT1.V());
3085 else if ( (beta<MyConfusionPrecision)&&(Abs(1-alpha)<MyConfusionPrecision) ) {//beta==0 alpha==1
3086 SP1.SetXYZ(PT2.X(),PT2.Y(),PT2.Z());
3087 SP1.SetUV2(PT2.U(),PT2.V());
3090 else if ( (Abs(beta-1)<MyConfusionPrecision)&&(Abs(1-alpha)<MyConfusionPrecision) ) {//beta==1 alpha==1
3091 SP1.SetXYZ(PT3.X(),PT3.Y(),PT3.Z());
3092 SP1.SetUV2(PT3.U(),PT3.V());
3095 else if (beta<MyConfusionPrecision) { //beta==0
3096 SP1.SetEdge2(Tri2.GetEdgeNumber(1));
3097 if(Tri2.GetEdgeOrientation(1)>0)
3098 SP1.SetLambda2(alpha);
3099 else SP1.SetLambda2(1.0-alpha);
3101 else if (Abs(beta-alpha)<MyConfusionPrecision) {//beta==alpha
3102 SP1.SetEdge2(Tri2.GetEdgeNumber(3));
3103 if(Tri2.GetEdgeOrientation(3)>0)
3104 SP1.SetLambda2(1.0-alpha);
3105 else SP1.SetLambda2(alpha);
3107 else if (Abs(alpha-1)<MyConfusionPrecision) {//alpha==1
3108 SP1.SetEdge2(Tri2.GetEdgeNumber(2));
3109 if(Tri2.GetEdgeOrientation(2)>0)
3110 SP1.SetLambda2(alpha);
3111 else SP1.SetLambda2(1.0-alpha);
3123 //=======================================================================
3124 //function : TriangleComparePSP
3125 //purpose : The same as TriangleCompare, plus compute the
3126 // StartPoints without chaining them.
3127 //=======================================================================
3128 Standard_Integer IntPolyh_MaillageAffinage::TriangleComparePSP ()
3130 Standard_Integer CpteurTab=0;
3131 Standard_Integer CpteurTabSP=0;
3132 Standard_Real CoupleAngle=-2.0;
3133 const Standard_Integer FinTT1 = TTriangles1.NbItems();
3134 const Standard_Integer FinTT2 = TTriangles2.NbItems();
3136 for(Standard_Integer i_S1=0; i_S1<FinTT1; i_S1++) {
3137 IntPolyh_Triangle &Triangle1 = TTriangles1[i_S1];
3138 if ((Triangle1.IndiceIntersectionPossible() == 0) ||
3139 (Triangle1.GetFleche() < 0.))
3141 for(Standard_Integer i_S2=0; i_S2<FinTT2; i_S2++){
3142 IntPolyh_Triangle &Triangle2 = TTriangles2[i_S2];
3143 if ((Triangle2.IndiceIntersectionPossible() != 0) &&
3144 (Triangle2.GetFleche() >= 0.)) {
3145 IntPolyh_StartPoint SP1, SP2;
3146 //If a triangle is dead or not in BSB, comparison is not possible
3148 Standard_Integer iDeg1, iDeg2, iDeg3, iDeg;
3150 const IntPolyh_Point& P1=TPoints1[Triangle1.FirstPoint()];
3151 const IntPolyh_Point& P2=TPoints1[Triangle1.SecondPoint()];
3152 const IntPolyh_Point& P3=TPoints1[Triangle1.ThirdPoint()];
3153 iDeg1=(P1.Degenerated()) ? 1 : 0;
3154 iDeg2=(P2.Degenerated()) ? 1 : 0;
3155 iDeg3=(P3.Degenerated()) ? 1 : 0;
3156 iDeg=iDeg1+iDeg2+iDeg3;
3161 const IntPolyh_Point& Q1=TPoints2[Triangle2.FirstPoint()];
3162 const IntPolyh_Point& Q2=TPoints2[Triangle2.SecondPoint()];
3163 const IntPolyh_Point& Q3=TPoints2[Triangle2.ThirdPoint()];
3164 iDeg1=(Q1.Degenerated()) ? 1 : 0;
3165 iDeg2=(Q2.Degenerated()) ? 1 : 0;
3166 iDeg3=(Q3.Degenerated()) ? 1 : 0;
3167 iDeg=iDeg1+iDeg2+iDeg3;
3172 if (TriContact(P1, P2, P3, Q1, Q2, Q3, CoupleAngle)) {
3173 Triangle1.SetIndiceIntersection(1);//The triangle is cut by another
3174 Triangle2.SetIndiceIntersection(1);
3176 Standard_Integer NbPoints;
3177 NbPoints=StartingPointsResearch(i_S1,i_S2,SP1, SP2);
3183 if ( (NbPoints>0)&&(NbPoints<3) ) {
3184 SP1.SetCoupleValue(i_S1,i_S2);
3185 TStartPoints[CpteurTabSP]=SP1;
3192 SP2.SetCoupleValue(i_S1,i_S2);
3193 TStartPoints[CpteurTabSP]=SP2;
3207 return(CpteurTabSP);
3209 //=======================================================================
3210 //function : TriangleCompare
3211 //purpose : Analyze each couple of triangles from the two --
3212 // array of triangles, to see if they are in
3213 // contact, and compute the incidence. Then put
3214 // couples in contact in the array of couples
3215 //=======================================================================
3216 Standard_Integer IntPolyh_MaillageAffinage::TriangleCompare ()
3218 Standard_Integer CpteurTab=0;
3220 const Standard_Integer FinTT1 = TTriangles1.NbItems();
3221 const Standard_Integer FinTT2 = TTriangles2.NbItems();
3223 Standard_Integer TTClimit = 200;
3224 Standard_Integer NbTTC = FinTT1 * FinTT2 / 10;
3225 if (NbTTC < TTClimit)
3227 TTrianglesContacts.Init(NbTTC);
3229 //TTrianglesContacts.Init(FinTT1 * FinTT2 / 10);
3231 Standard_Real CoupleAngle=-2.0;
3232 for(Standard_Integer i_S1=0; i_S1<FinTT1; i_S1++) {
3233 IntPolyh_Triangle &Triangle1 = TTriangles1[i_S1];
3234 if ((Triangle1.IndiceIntersectionPossible() == 0) ||
3235 (Triangle1.GetFleche() < 0.))
3237 for(Standard_Integer i_S2=0; i_S2<FinTT2; i_S2++){
3238 IntPolyh_Triangle &Triangle2 = TTriangles2[i_S2];
3239 if ((Triangle2.IndiceIntersectionPossible() != 0) &&
3240 (Triangle2.GetFleche() >= 0.)) {
3241 //If a triangle is dead or not in BSB, comparison is not possible
3242 Standard_Integer iDeg1, iDeg2, iDeg3, iDeg;
3244 const IntPolyh_Point& P1=TPoints1[Triangle1.FirstPoint()];
3245 const IntPolyh_Point& P2=TPoints1[Triangle1.SecondPoint()];
3246 const IntPolyh_Point& P3=TPoints1[Triangle1.ThirdPoint()];
3247 iDeg1=(P1.Degenerated()) ? 1 : 0;
3248 iDeg2=(P2.Degenerated()) ? 1 : 0;
3249 iDeg3=(P3.Degenerated()) ? 1 : 0;
3250 iDeg=iDeg1+iDeg2+iDeg3;
3255 const IntPolyh_Point& Q1=TPoints2[Triangle2.FirstPoint()];
3256 const IntPolyh_Point& Q2=TPoints2[Triangle2.SecondPoint()];
3257 const IntPolyh_Point& Q3=TPoints2[Triangle2.ThirdPoint()];
3258 iDeg1=(Q1.Degenerated()) ? 1 : 0;
3259 iDeg2=(Q2.Degenerated()) ? 1 : 0;
3260 iDeg3=(Q3.Degenerated()) ? 1 : 0;
3261 iDeg=iDeg1+iDeg2+iDeg3;
3266 if (TriContact(P1, P2, P3, Q1, Q2, Q3, CoupleAngle)) {
3267 if (CpteurTab >= NbTTC)
3269 TTrianglesContacts.SetNbItems(CpteurTab);
3272 TTrianglesContacts[CpteurTab].SetCoupleValue(i_S1, i_S2);
3273 TTrianglesContacts[CpteurTab].SetAngleValue(CoupleAngle);
3274 //test TTrianglesContacts[CpteurTab].Dump(CpteurTab);
3276 Triangle1.SetIndiceIntersection(1);//The triangle is cut by another
3277 Triangle2.SetIndiceIntersection(1);
3283 TTrianglesContacts.SetNbItems(CpteurTab);
3288 //=======================================================================
3289 //function : StartPointsCalcul
3290 //purpose : From the array of couples compute all the start
3291 // points and display them on the screen
3292 //=======================================================================
3293 void IntPolyh_MaillageAffinage::StartPointsCalcul() const
3295 const Standard_Integer FinTTC = TTrianglesContacts.NbItems();
3296 // printf("StartPointsCalcul() from IntPolyh_MaillageAffinage.cxx : StartPoints:\n");
3297 for(Standard_Integer ii=0; ii<FinTTC; ii++) {
3298 IntPolyh_StartPoint SP1,SP2;
3299 Standard_Integer T1,T2;
3300 T1=TTrianglesContacts[ii].FirstValue();
3301 T2=TTrianglesContacts[ii].SecondValue();
3302 StartingPointsResearch(T1,T2,SP1,SP2);
3303 if ( (SP1.E1()!=-1)&&(SP1.E2()!=-1) ) SP1.Dump(ii);
3304 if ( (SP2.E1()!=-1)&&(SP2.E2()!=-1) ) SP2.Dump(ii);
3307 //=======================================================================
3308 //function : CheckCoupleAndGetAngle
3310 //=======================================================================
3311 Standard_Integer CheckCoupleAndGetAngle(const Standard_Integer T1,
3312 const Standard_Integer T2,
3313 Standard_Real& Angle,
3314 IntPolyh_ArrayOfCouples &TTrianglesContacts)
3316 Standard_Integer Test=0;
3317 const Standard_Integer FinTTC = TTrianglesContacts.NbItems();
3318 for (Standard_Integer oioi=0; oioi<FinTTC; oioi++) {
3319 IntPolyh_Couple TestCouple = TTrianglesContacts[oioi];
3320 if ( (TestCouple.FirstValue()==T1)&&(TestCouple.AnalyseFlagValue()!=1) ) {
3321 if (TestCouple.SecondValue()==T2) {
3323 TTrianglesContacts[oioi].SetAnalyseFlag(1);
3324 Angle=TTrianglesContacts[oioi].AngleValue();
3331 //=======================================================================
3332 //function : CheckCoupleAndGetAngle2
3334 //=======================================================================
3335 Standard_Integer CheckCoupleAndGetAngle2(const Standard_Integer T1,
3336 const Standard_Integer T2,
3337 const Standard_Integer T11,
3338 const Standard_Integer T22,
3339 Standard_Integer &CT11,
3340 Standard_Integer &CT22,
3341 Standard_Real & Angle,
3342 IntPolyh_ArrayOfCouples &TTrianglesContacts)
3344 ///couple T1 T2 is found in the list
3345 ///T11 and T22 are two other triangles implied in the contact edge edge
3346 /// CT11 couple( T1,T22) and CT22 couple (T2,T11)
3347 /// these couples will be marked if there is a start point
3348 Standard_Integer Test1=0;
3349 Standard_Integer Test2=0;
3350 Standard_Integer Test3=0;
3351 const Standard_Integer FinTTC = TTrianglesContacts.NbItems();
3352 for (Standard_Integer oioi=0; oioi<FinTTC; oioi++) {
3353 IntPolyh_Couple TestCouple = TTrianglesContacts[oioi];
3354 if( (Test1==0)||(Test2==0)||(Test3==0) ) {
3355 if ( (TestCouple.FirstValue()==T1)&&(TestCouple.AnalyseFlagValue()!=1) ) {
3356 if (TestCouple.SecondValue()==T2) {
3358 TTrianglesContacts[oioi].SetAnalyseFlag(1);
3359 Angle=TTrianglesContacts[oioi].AngleValue();
3361 else if (TestCouple.SecondValue()==T22) {
3364 Angle=TTrianglesContacts[oioi].AngleValue();
3367 else if( (TestCouple.FirstValue()==T11)&&(TestCouple.AnalyseFlagValue()!=1) ) {
3368 if (TestCouple.SecondValue()==T2) {
3371 Angle=TTrianglesContacts[oioi].AngleValue();
3380 //=======================================================================
3381 //function : CheckNextStartPoint
3382 //purpose : it is checked if the point is not a top
3383 // then it is stored in one or several valid arrays with
3384 // the proper list number
3385 //=======================================================================
3386 Standard_Integer CheckNextStartPoint(IntPolyh_SectionLine & SectionLine,
3387 IntPolyh_ArrayOfTangentZones & TTangentZones,
3388 IntPolyh_StartPoint & SP,
3389 const Standard_Boolean Prepend)//=Standard_False)
3391 Standard_Integer Test=1;
3392 if( (SP.E1()==-1)||(SP.E2()==-1) ) {
3393 //The tops of triangle are analyzed
3394 //It is checked if they are not in the array TTangentZones
3395 Standard_Integer FinTTZ=TTangentZones.NbItems();
3396 for(Standard_Integer uiui=0; uiui<FinTTZ; uiui++) {
3397 IntPolyh_StartPoint TestSP=TTangentZones[uiui];
3398 if ( (Abs(SP.U1()-TestSP.U1())<MyConfusionPrecision)
3399 &&(Abs(SP.V1()-TestSP.V1())<MyConfusionPrecision) ) {
3400 if ( (Abs(SP.U2()-TestSP.U2())<MyConfusionPrecision)
3401 &&(Abs(SP.V2()-TestSP.V2())<MyConfusionPrecision) ) {
3402 Test=0;//SP is already in the list of tops
3407 if (Test) {//the top does not belong to the list of TangentZones
3408 SP.SetChainList(-1);
3409 TTangentZones[FinTTZ]=SP;
3410 TTangentZones.IncrementNbItems();
3411 Test=0;//the examined point is a top
3416 SectionLine.Prepend(SP);
3418 SectionLine[SectionLine.NbStartPoints()]=SP;
3419 SectionLine.IncrementNbStartPoints();
3423 //if the point is not a top Test=1
3424 //The chain is continued
3427 //=======================================================================
3428 //function : StartPointsChain
3429 //purpose : Loop on the array of couples. Compute StartPoints.
3430 // Try to chain the StartPoints into SectionLines or
3431 // put the point in the ArrayOfTangentZones if
3432 // chaining it, is not possible.
3433 //=======================================================================
3434 Standard_Integer IntPolyh_MaillageAffinage::StartPointsChain
3435 (IntPolyh_ArrayOfSectionLines& TSectionLines,
3436 IntPolyh_ArrayOfTangentZones& TTangentZones)
3438 //Loop on the array of couples filled in the function COMPARE()
3439 const Standard_Integer FinTTC = TTrianglesContacts.NbItems();
3441 //Array of tops of triangles
3442 for(Standard_Integer IndexA=0; IndexA<FinTTC; IndexA++) {
3443 //First couple of triangles.
3444 //It is checked if the couple of triangles has not been already examined.
3445 if(TTrianglesContacts[IndexA].AnalyseFlagValue()!=1) {
3447 Standard_Integer SectionLineIndex=TSectionLines.NbItems();
3448 // fill last section line if still empty (eap)
3449 if (SectionLineIndex > 0
3451 TSectionLines[SectionLineIndex-1].NbStartPoints() == 0)
3452 SectionLineIndex -= 1;
3454 TSectionLines.IncrementNbItems();
3456 IntPolyh_SectionLine & MySectionLine=TSectionLines[SectionLineIndex];
3457 if (MySectionLine.GetN() == 0) // eap
3458 MySectionLine.Init(10000);//Initialisation of array of StartPoint
3460 Standard_Integer NbPoints=-1;
3461 Standard_Integer T1I, T2I;
3462 T1I = TTrianglesContacts[IndexA].FirstValue();
3463 T2I = TTrianglesContacts[IndexA].SecondValue();
3465 // Start points for the current couple are found
3466 IntPolyh_StartPoint SP1, SP2;
3467 NbPoints=StartingPointsResearch2(T1I,T2I,SP1, SP2);//first calculation
3468 TTrianglesContacts[IndexA].SetAnalyseFlag(1);//the couple is marked
3470 if(NbPoints==1) {// particular case top/triangle or edge/edge
3471 //the start point is input in the array
3472 SP1.SetChainList(SectionLineIndex);
3473 SP1.SetAngle(TTrianglesContacts[IndexA].AngleValue());
3474 //it is checked if the point is not atop of the triangle
3475 if(CheckNextStartPoint(MySectionLine,TTangentZones,SP1)) {
3476 IntPolyh_StartPoint SPNext1;
3477 Standard_Integer TestSP1=0;
3480 IntPolyh_StartPoint SP11;//=SP1;
3481 if(SP1.E1()>=0) { //&&(SP1.E2()!=-1) already tested if the point is not a top
3482 Standard_Integer NextTriangle1=0;
3483 if (TEdges1[SP1.E1()].FirstTriangle()!=T1I) NextTriangle1=TEdges1[SP1.E1()].FirstTriangle();
3484 else NextTriangle1=TEdges1[SP1.E1()].SecondTriangle();
3486 Standard_Real Angle=-2.0;
3487 if (CheckCoupleAndGetAngle(NextTriangle1,T2I,Angle,TTrianglesContacts)) {
3488 //it is checked if the couple exists and is marked
3489 Standard_Integer NbPoints11=0;
3490 NbPoints11=NextStartingPointsResearch2(NextTriangle1,T2I,SP1,SP11);
3491 if (NbPoints11==1) {
3492 SP11.SetChainList(SectionLineIndex);
3493 SP11.SetAngle(Angle);
3495 if(CheckNextStartPoint(MySectionLine,TTangentZones,SP11)) {
3496 Standard_Integer EndChainList=1;
3497 while (EndChainList!=0) {
3498 TestSP1=GetNextChainStartPoint(SP11,SPNext1,MySectionLine,TTangentZones);
3500 SPNext1.SetChainList(SectionLineIndex);
3501 if(CheckNextStartPoint(MySectionLine,TTangentZones,SPNext1))
3503 else EndChainList=0;
3505 else EndChainList=0; //There is no next point
3511 if(NbPoints11>1) {//The point is input in the array TTangentZones
3512 TTangentZones[TTangentZones.NbItems()]=SP11;//default list number = -1
3513 TTangentZones.IncrementNbItems();
3521 else if (SP1.E2()<0){
3524 //chain of the other side
3525 IntPolyh_StartPoint SP12;//=SP1;
3526 if (SP1.E2()>=0) { //&&(SP1.E1()!=-1) already tested
3527 Standard_Integer NextTriangle2;
3528 if (TEdges2[SP1.E2()].FirstTriangle()!=T2I) NextTriangle2=TEdges2[SP1.E2()].FirstTriangle();
3529 else NextTriangle2=TEdges2[SP1.E2()].SecondTriangle();
3531 Standard_Real Angle=-2.0;
3532 if(CheckCoupleAndGetAngle(T1I,NextTriangle2,Angle,TTrianglesContacts)) {
3533 Standard_Integer NbPoints12=0;
3534 NbPoints12=NextStartingPointsResearch2(T1I,NextTriangle2,SP1, SP12);
3535 if (NbPoints12==1) {
3537 SP12.SetChainList(SectionLineIndex);
3538 SP12.SetAngle(Angle);
3539 Standard_Boolean Prepend = Standard_True; // eap
3541 if(CheckNextStartPoint(MySectionLine,TTangentZones,SP12, Prepend)) {
3542 Standard_Integer EndChainList=1;
3543 while (EndChainList!=0) {
3544 TestSP1=GetNextChainStartPoint(SP12,SPNext1,
3545 MySectionLine,TTangentZones,
3548 SPNext1.SetChainList(SectionLineIndex);
3549 if(CheckNextStartPoint(MySectionLine,TTangentZones,SPNext1,Prepend))
3551 else EndChainList=0;
3553 else EndChainList=0; //there is no next point
3558 if(NbPoints12>1) {//The points are input in the array TTangentZones
3559 TTangentZones[TTangentZones.NbItems()]=SP12;//default list number = -1
3560 TTangentZones.IncrementNbItems();
3569 else if(SP1.E1()<0){
3574 else if(NbPoints==2) {
3575 //the start points are input in the array
3576 IntPolyh_StartPoint SPNext2;
3577 Standard_Integer TestSP2=0;
3578 Standard_Integer EndChainList=1;
3580 SP1.SetChainList(SectionLineIndex);
3581 SP1.SetAngle(TTrianglesContacts[IndexA].AngleValue());
3582 if(CheckNextStartPoint(MySectionLine,TTangentZones,SP1)) {
3585 while (EndChainList!=0) {
3586 TestSP2=GetNextChainStartPoint(SP1,SPNext2,MySectionLine,TTangentZones);
3588 SPNext2.SetChainList(SectionLineIndex);
3589 if(CheckNextStartPoint(MySectionLine,TTangentZones,SPNext2))
3591 else EndChainList=0;
3593 else EndChainList=0; //there is no next point
3597 SP2.SetChainList(SectionLineIndex);
3598 SP2.SetAngle(TTrianglesContacts[IndexA].AngleValue());
3599 Standard_Boolean Prepend = Standard_True; // eap
3601 if(CheckNextStartPoint(MySectionLine,TTangentZones,SP2,Prepend)) {
3603 //chain of the other side
3605 while (EndChainList!=0) {
3606 TestSP2=GetNextChainStartPoint(SP2,SPNext2,
3607 MySectionLine,TTangentZones,
3610 SPNext2.SetChainList(SectionLineIndex);
3611 if(CheckNextStartPoint(MySectionLine,TTangentZones,SPNext2,Prepend))
3613 else EndChainList=0;
3615 else EndChainList=0; //there is no next point
3620 else if( (NbPoints>2)&&(NbPoints<7) ) {
3621 //More than two start points
3622 //the start point is input in the table
3623 SP1.SetChainList(SectionLineIndex);
3624 CheckNextStartPoint(MySectionLine,TTangentZones,SP1);
3635 //=======================================================================
3636 //function : GetNextChainStartPoint
3637 //purpose : Mainly used by StartPointsChain(), this function
3638 // try to compute the next StartPoint.
3639 // GetNextChainStartPoint is used only if it is known that there are 2 contact points
3640 //=======================================================================
3641 Standard_Integer IntPolyh_MaillageAffinage::GetNextChainStartPoint
3642 (const IntPolyh_StartPoint & SP,
3643 IntPolyh_StartPoint & SPNext,
3644 IntPolyh_SectionLine & MySectionLine,
3645 IntPolyh_ArrayOfTangentZones & TTangentZones,
3646 const Standard_Boolean Prepend)
3648 Standard_Integer NbPoints=0;
3649 if( (SP.E1()>=0)&&(SP.E2()==-2) ) {
3650 //case if the point is on edge of T1
3651 Standard_Integer NextTriangle1;
3652 if (TEdges1[SP.E1()].FirstTriangle()!=SP.T1()) NextTriangle1=TEdges1[SP.E1()].FirstTriangle();
3654 NextTriangle1=TEdges1[SP.E1()].SecondTriangle();
3655 //If is checked if two triangles intersect
3656 Standard_Real Angle= -2.0;
3657 if (CheckCoupleAndGetAngle(NextTriangle1,SP.T2(),Angle,TTrianglesContacts)) {
3658 NbPoints=NextStartingPointsResearch2(NextTriangle1,SP.T2(),SP,SPNext);
3661 CheckNextStartPoint(MySectionLine,TTangentZones,SPNext,Prepend);
3668 SPNext.SetAngle(Angle);
3670 else NbPoints=0;//this couple does not intersect
3672 else if( (SP.E1()==-2)&&(SP.E2()>=0) ) {
3673 //case if the point is on edge of T2
3674 Standard_Integer NextTriangle2;
3675 if (TEdges2[SP.E2()].FirstTriangle()!=SP.T2()) NextTriangle2=TEdges2[SP.E2()].FirstTriangle();
3677 NextTriangle2=TEdges2[SP.E2()].SecondTriangle();
3678 Standard_Real Angle= -2.0;
3679 if (CheckCoupleAndGetAngle(SP.T1(),NextTriangle2,Angle,TTrianglesContacts)) {
3680 NbPoints=NextStartingPointsResearch2(SP.T1(),NextTriangle2,SP,SPNext);
3683 CheckNextStartPoint(MySectionLine,TTangentZones,SPNext,Prepend);
3690 SPNext.SetAngle(Angle);
3694 else if( (SP.E1()==-2)&&(SP.E2()==-2) ) {
3695 ///no edge is touched or cut
3699 else if( (SP.E1()>=0)&&(SP.E2()>=0) ) {
3700 ///the point is located on two edges
3701 Standard_Integer NextTriangle1;
3702 Standard_Integer CpleT11=-1;
3703 Standard_Integer CpleT22=-1;
3704 if (TEdges1[SP.E1()].FirstTriangle()!=SP.T1()) NextTriangle1=TEdges1[SP.E1()].FirstTriangle();
3706 NextTriangle1=TEdges1[SP.E1()].SecondTriangle();
3707 Standard_Integer NextTriangle2;
3708 if (TEdges2[SP.E2()].FirstTriangle()!=SP.T2()) NextTriangle2=TEdges2[SP.E2()].FirstTriangle();
3710 NextTriangle2=TEdges2[SP.E2()].SecondTriangle();
3711 Standard_Real Angle= -2.0;
3712 if (CheckCoupleAndGetAngle2(NextTriangle1,NextTriangle2,
3713 SP.T1(),SP.T2(),CpleT11,CpleT22,
3714 Angle,TTrianglesContacts)) {
3715 NbPoints=NextStartingPointsResearch2(NextTriangle1,NextTriangle2,SP,SPNext);
3718 ///The new point is checked
3719 if(CheckNextStartPoint(MySectionLine,TTangentZones,SPNext,Prepend)>0) {
3728 SPNext.SetAngle(Angle);
3729 //The couples (Ti,Tj) (Ti',Tj') are marked
3730 if (CpleT11>=0) TTrianglesContacts[CpleT11].SetAnalyseFlag(1);
3734 if (CpleT22>=0) TTrianglesContacts[CpleT22].SetAnalyseFlag(1);
3742 else if( (SP.E1()==-1)||(SP.E2()==-1) ) {
3743 ///the points are tops of triangle
3744 ///the point is atored in an intermediary array
3748 //=======================================================================
3749 //function : GetArrayOfPoints
3751 //=======================================================================
3752 const IntPolyh_ArrayOfPoints& IntPolyh_MaillageAffinage::GetArrayOfPoints
3753 (const Standard_Integer SurfID)const
3759 //=======================================================================
3760 //function : GetArrayOfEdges
3762 //=======================================================================
3763 const IntPolyh_ArrayOfEdges& IntPolyh_MaillageAffinage::GetArrayOfEdges
3764 (const Standard_Integer SurfID)const
3770 //=======================================================================
3771 //function : GetArrayOfTriangles
3773 //=======================================================================
3774 const IntPolyh_ArrayOfTriangles&
3775 IntPolyh_MaillageAffinage::GetArrayOfTriangles
3776 (const Standard_Integer SurfID)const{
3778 return(TTriangles1);
3779 return(TTriangles2);
3782 //=======================================================================
3785 //=======================================================================
3786 Bnd_Box IntPolyh_MaillageAffinage::GetBox(const Standard_Integer SurfID) const
3793 //=======================================================================
3794 //function : GetBoxDraw
3796 //=======================================================================
3797 void IntPolyh_MaillageAffinage::GetBoxDraw(const Standard_Integer SurfID)const
3799 Standard_Real x0,y0,z0,x1,y1,z1;
3801 MyBox1.Get(x0,y0,z0,x1,y1,z1);
3804 MyBox2.Get(x0,y0,z0,x1,y1,z1);
3807 //=======================================================================
3808 //function : GetArrayOfCouples
3810 //=======================================================================
3811 IntPolyh_ArrayOfCouples &IntPolyh_MaillageAffinage::GetArrayOfCouples()
3813 return TTrianglesContacts;
3815 //=======================================================================
3816 //function : SetEnlargeZone
3818 //=======================================================================
3819 void IntPolyh_MaillageAffinage::SetEnlargeZone(Standard_Boolean& EnlargeZone)
3821 myEnlargeZone = EnlargeZone;
3823 //=======================================================================
3824 //function : GetEnlargeZone
3826 //=======================================================================
3827 Standard_Boolean IntPolyh_MaillageAffinage::GetEnlargeZone() const
3829 return myEnlargeZone;
3832 //=======================================================================
3833 //function : GetMinDeflection
3835 //=======================================================================
3836 Standard_Real IntPolyh_MaillageAffinage::GetMinDeflection(const Standard_Integer SurfID) const
3838 return (SurfID==1)? FlecheMin1:FlecheMin2;
3841 //=======================================================================
3842 //function : GetMaxDeflection
3844 //=======================================================================
3845 Standard_Real IntPolyh_MaillageAffinage::GetMaxDeflection(const Standard_Integer SurfID) const
3847 return (SurfID==1)? FlecheMax1:FlecheMax2;
3850 //=======================================================================
3851 //function : DegeneratedIndex
3853 //=======================================================================
3854 void DegeneratedIndex(const TColStd_Array1OfReal& aXpars,
3855 const Standard_Integer aNbX,
3856 const Handle(Adaptor3d_HSurface)& aS,
3857 const Standard_Integer aIsoDirection,
3858 Standard_Integer& aI1,
3859 Standard_Integer& aI2)
3862 Standard_Boolean bDegX1, bDegX2;
3863 Standard_Real aDegX1, aDegX2, aTol2, aX;
3867 aTol2=MyTolerance*MyTolerance;
3869 if (aIsoDirection==1){ // V=const
3870 bDegX1=IsDegenerated(aS, 1, aTol2, aDegX1);
3871 bDegX2=IsDegenerated(aS, 2, aTol2, aDegX2);
3873 else if (aIsoDirection==2){ // U=const
3874 bDegX1=IsDegenerated(aS, 3, aTol2, aDegX1);
3875 bDegX2=IsDegenerated(aS, 4, aTol2, aDegX2);
3881 if (!(bDegX1 || bDegX2)) {
3885 for(i=1; i<=aNbX; ++i) {
3888 if (fabs(aX-aDegX1) < MyTolerance) {
3893 if (fabs(aX-aDegX2) < MyTolerance) {
3899 //=======================================================================
3900 //function : IsDegenerated
3902 //=======================================================================
3903 Standard_Boolean IsDegenerated(const Handle(Adaptor3d_HSurface)& aS,
3904 const Standard_Integer aIndex,
3905 const Standard_Real aTol2,
3906 Standard_Real& aDegX)
3908 Standard_Boolean bRet;
3909 Standard_Integer i, aNbP;
3910 Standard_Real aU, dU, aU1, aU2, aV, dV, aV1, aV2, aD2;
3913 bRet=Standard_False;
3917 aU1=aS->FirstUParameter();
3918 aU2=aS->LastUParameter();
3919 aV1=aS->FirstVParameter();
3920 aV2=aS->LastVParameter();
3922 if (aIndex<3) { // V=const
3927 dU=(aU2-aU1)/(aNbP-1);
3929 aP1=aS->Value(aU, aV);
3930 for (i=1; i<aNbP; ++i) {
3935 aP2=aS->Value(aU, aV);
3936 aD2=aP1.SquareDistance(aP2);
3950 dV=(aV2-aV1)/(aNbP-1);
3952 aP1=aS->Value(aU, aV);
3953 for (i=1; i<aNbP; ++i) {
3958 aP2=aS->Value(aU, aV);
3959 aD2=aP1.SquareDistance(aP2);
3971 //=======================================================================
3972 //function : EnlargeZone
3974 //=======================================================================
3975 void EnlargeZone(const Handle(Adaptor3d_HSurface)& MaSurface,
3981 if(MaSurface->GetType() == GeomAbs_BSplineSurface ||
3982 MaSurface->GetType() == GeomAbs_BezierSurface) {
3983 if((!MaSurface->IsUClosed() && !MaSurface->IsUPeriodic()) &&
3984 (Abs(u0) < 1.e+100 && Abs(u1) < 1.e+100) ) {
3985 Standard_Real delta_u = 0.01*Abs(u1 - u0);
3989 if((!MaSurface->IsVClosed() && !MaSurface->IsVPeriodic()) &&
3990 (Abs(v0) < 1.e+100 && Abs(v1) < 1.e+100) ) {
3991 Standard_Real delta_v = 0.01*Abs(v1 - v0);