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
9 // under the terms of the GNU Lesser General Public 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
21 // modified by Edward AGAPOV (eap) Thu Feb 14 2002 (occ139)
22 // - make Section Line parts rightly connected (prepend 2nd part to the 1st)
23 // - TriangleShape() for debugging purpose
25 // Modified by skv - Thu Sep 25 17:42:42 2003 OCC567
26 // modified by ofv Thu Apr 8 14:58:13 2004 fip
29 #include <IntPolyh_MaillageAffinage.ixx>
31 #include <Precision.hxx>
35 #include <TColStd_ListIteratorOfListOfInteger.hxx>
37 #include <Bnd_BoundSortBox.hxx>
38 #include <Bnd_HArray1OfBox.hxx>
40 #include <IntCurveSurface_ThePolyhedronOfHInter.hxx>
42 #include <IntPolyh_ArrayOfCouples.hxx>
43 #include <IntPolyh_Edge.hxx>
44 #include <IntPolyh_Couple.hxx>
46 static Standard_Real MyTolerance=10.0e-7;
47 static Standard_Real MyConfusionPrecision=10.0e-12;
48 static Standard_Real SquareMyConfusionPrecision=10.0e-24;
51 inline Standard_Real maxSR(const Standard_Real a,
52 const Standard_Real b,
53 const Standard_Real c);
56 inline Standard_Real minSR(const Standard_Real a,
57 const Standard_Real b,
58 const Standard_Real c);
60 Standard_Integer project6(const IntPolyh_Point &ax,
61 const IntPolyh_Point &p1,
62 const IntPolyh_Point &p2,
63 const IntPolyh_Point &p3,
64 const IntPolyh_Point &q1,
65 const IntPolyh_Point &q2,
66 const IntPolyh_Point &q3);
68 void TestNbPoints(const Standard_Integer ,
69 Standard_Integer &NbPoints,
70 Standard_Integer &NbPointsTotal,
71 const IntPolyh_StartPoint &Pt1,
72 const IntPolyh_StartPoint &Pt2,
73 IntPolyh_StartPoint &SP1,
74 IntPolyh_StartPoint &SP2);
76 void CalculPtsInterTriEdgeCoplanaires(const Standard_Integer TriSurfID,
77 const IntPolyh_Point &NormaleTri,
78 const IntPolyh_Point &PE1,
79 const IntPolyh_Point &PE2,
80 const IntPolyh_Point &Edge,
81 const IntPolyh_Point &PT1,
82 const IntPolyh_Point &PT2,
83 const IntPolyh_Point &Cote,
84 const Standard_Integer CoteIndex,
85 IntPolyh_StartPoint &SP1,
86 IntPolyh_StartPoint &SP2,
87 Standard_Integer &NbPoints);
89 void CalculPtsInterTriEdgeCoplanaires2(const Standard_Integer TriSurfID,
90 const IntPolyh_Point &NormaleTri,
91 const IntPolyh_Triangle &Tri1,
92 const IntPolyh_Triangle &Tri2,
93 const IntPolyh_Point &PE1,
94 const IntPolyh_Point &PE2,
95 const IntPolyh_Point &Edge,
96 const Standard_Integer EdgeIndex,
97 const IntPolyh_Point &PT1,
98 const IntPolyh_Point &PT2,
99 const IntPolyh_Point &Cote,
100 const Standard_Integer CoteIndex,
101 IntPolyh_StartPoint &SP1,
102 IntPolyh_StartPoint &SP2,
103 Standard_Integer &NbPoints);
105 Standard_Integer CheckCoupleAndGetAngle(const Standard_Integer T1,
106 const Standard_Integer T2,
107 Standard_Real& Angle,
108 IntPolyh_ArrayOfCouples &TTrianglesContacts);
110 Standard_Integer CheckCoupleAndGetAngle2(const Standard_Integer T1,
111 const Standard_Integer T2,
112 const Standard_Integer T11,
113 const Standard_Integer T22,
114 Standard_Integer &CT11,
115 Standard_Integer &CT22,
116 Standard_Real & Angle,
117 IntPolyh_ArrayOfCouples &TTrianglesContacts);
119 Standard_Integer CheckNextStartPoint(IntPolyh_SectionLine & SectionLine,
120 IntPolyh_ArrayOfTangentZones & TTangentZones,
121 IntPolyh_StartPoint & SP,
122 const Standard_Boolean Prepend=Standard_False);
125 Standard_Boolean IsDegenerated(const Handle(Adaptor3d_HSurface)& aS,
126 const Standard_Integer aIndex,
127 const Standard_Real aTol2,
128 Standard_Real& aDegX);
130 void DegeneratedIndex(const TColStd_Array1OfReal& Xpars,
131 const Standard_Integer aNbX,
132 const Handle(Adaptor3d_HSurface)& aS,
133 const Standard_Integer aIsoDirection,
134 Standard_Integer& aI1,
135 Standard_Integer& aI2);
137 void EnlargeZone(const Handle(Adaptor3d_HSurface)& MaSurface,
143 //=======================================================================
144 //function : IntPolyh_MaillageAffinage
146 //=======================================================================
147 IntPolyh_MaillageAffinage::IntPolyh_MaillageAffinage
148 (const Handle(Adaptor3d_HSurface)& Surface1,
149 const Handle(Adaptor3d_HSurface)& Surface2,
150 const Standard_Integer )
152 MaSurface1(Surface1),
153 MaSurface2(Surface2),
164 myEnlargeZone(Standard_False)
167 //=======================================================================
168 //function : IntPolyh_MaillageAffinage
170 //=======================================================================
171 IntPolyh_MaillageAffinage::IntPolyh_MaillageAffinage
172 (const Handle(Adaptor3d_HSurface)& Surface1,
173 const Standard_Integer NbSU1,
174 const Standard_Integer NbSV1,
175 const Handle(Adaptor3d_HSurface)& Surface2,
176 const Standard_Integer NbSU2,
177 const Standard_Integer NbSV2,
178 const Standard_Integer )
180 MaSurface1(Surface1),
181 MaSurface2(Surface2),
192 myEnlargeZone(Standard_False)
195 //=======================================================================
196 //function : FillArrayOfPnt
197 //purpose : Compute points on one surface and fill an array of points
198 //=======================================================================
199 void IntPolyh_MaillageAffinage::FillArrayOfPnt
200 (const Standard_Integer SurfID)
202 Standard_Integer NbSamplesU, NbSamplesV, i, aNbSamplesU1, aNbSamplesV1;
203 Standard_Real u0, u1, v0, v1, aU, aV, dU, dV;
205 const Handle(Adaptor3d_HSurface&) MaSurface=(SurfID==1)? MaSurface1 : MaSurface2;
206 NbSamplesU=(SurfID==1)? NbSamplesU1:NbSamplesU2;
207 NbSamplesV=(SurfID==1)? NbSamplesV1:NbSamplesV2;
209 u0 = (MaSurface)->FirstUParameter();
210 u1 = (MaSurface)->LastUParameter();
211 v0 = (MaSurface)->FirstVParameter();
212 v1 = (MaSurface)->LastVParameter();
215 EnlargeZone(MaSurface, u0, u1, v0, v1);
218 TColStd_Array1OfReal aUpars(1, NbSamplesU);
219 TColStd_Array1OfReal aVpars(1, NbSamplesV);
221 aNbSamplesU1=NbSamplesU-1;
222 aNbSamplesV1=NbSamplesV-1;
224 dU=(u1-u0)/Standard_Real(aNbSamplesU1);
225 dV=(v1-v0)/Standard_Real(aNbSamplesV1);
227 for (i=0; i<NbSamplesU; ++i) {
229 if (i==aNbSamplesU1) {
232 aUpars.SetValue(i+1, aU);
235 for (i=0; i<NbSamplesV; ++i) {
237 if (i==aNbSamplesV1) {
240 aVpars.SetValue(i+1, aV);
243 FillArrayOfPnt(SurfID, aUpars, aVpars);
245 //=======================================================================
246 //function : FillArrayOfPnt
247 //purpose : Compute points on one surface and fill an array of points
248 // FILL AN ARRAY OF POINTS
249 //=======================================================================
250 void IntPolyh_MaillageAffinage::FillArrayOfPnt
251 (const Standard_Integer SurfID,
252 const Standard_Boolean isShiftFwd)
254 Standard_Integer NbSamplesU, NbSamplesV, i, aNbSamplesU1, aNbSamplesV1;
255 Standard_Real u0, u1, v0, v1, aU, aV, dU, dV;
256 const Handle(Adaptor3d_HSurface)& MaSurface=(SurfID==1)? MaSurface1 : MaSurface2;
257 NbSamplesU=(SurfID==1)? NbSamplesU1:NbSamplesU2;
258 NbSamplesV=(SurfID==1)? NbSamplesV1:NbSamplesV2;
260 u0 = (MaSurface)->FirstUParameter();
261 u1 = (MaSurface)->LastUParameter();
262 v0 = (MaSurface)->FirstVParameter();
263 v1 = (MaSurface)->LastVParameter();
266 EnlargeZone(MaSurface, u0, u1, v0, v1);
269 TColStd_Array1OfReal aUpars(1, NbSamplesU);
270 TColStd_Array1OfReal aVpars(1, NbSamplesV);
272 aNbSamplesU1=NbSamplesU-1;
273 aNbSamplesV1=NbSamplesV-1;
275 dU=(u1-u0)/Standard_Real(aNbSamplesU1);
276 dV=(v1-v0)/Standard_Real(aNbSamplesV1);
278 for (i=0; i<NbSamplesU; ++i) {
280 if (i==aNbSamplesU1) {
283 aUpars.SetValue(i+1, aU);
286 for (i=0; i<NbSamplesV; ++i) {
288 if (i==aNbSamplesV1) {
291 aVpars.SetValue(i+1, aV);
294 FillArrayOfPnt(SurfID, isShiftFwd, aUpars, aVpars);
296 //=======================================================================
297 //function : FillArrayOfPnt
298 //purpose : Compute points on one surface and fill an array of points
299 //=======================================================================
300 void IntPolyh_MaillageAffinage::FillArrayOfPnt
301 (const Standard_Integer SurfID,
302 const TColStd_Array1OfReal& Upars,
303 const TColStd_Array1OfReal& Vpars)
305 Standard_Boolean bDegI, bDeg;
306 Standard_Integer aNbU, aNbV, iCnt, i, j;
307 Standard_Integer aID1, aID2, aJD1, aJD2;
308 Standard_Real aTol, aU, aV, aX, aY, aZ;
311 aNbU=(SurfID==1)? NbSamplesU1 : NbSamplesU2;
312 aNbV=(SurfID==1)? NbSamplesV1 : NbSamplesV2;
313 Bnd_Box& aBox = (SurfID==1) ? MyBox1 : MyBox2;
314 Handle(Adaptor3d_HSurface)& aS=(SurfID==1)? MaSurface1:MaSurface2;
315 IntPolyh_ArrayOfPoints &TPoints=(SurfID==1)? TPoints1:TPoints2;
321 DegeneratedIndex(Vpars, aNbV, aS, 1, aJD1, aJD2);
322 if (!(aJD1 || aJD2)) {
323 DegeneratedIndex(Upars, aNbU, aS, 2, aID1, aID2);
326 TPoints.Init(aNbU*aNbV);
328 for(i=1; i<=aNbU; ++i){
329 bDegI=(aID1==i || aID2==i);
331 for(j=1; j<=aNbV; ++j){
333 aP=aS->Value(aU, aV);
334 aP.Coord(aX, aY, aZ);
335 IntPolyh_Point& aIP=TPoints[iCnt];
336 aIP.Set(aX, aY, aZ, aU, aV);
338 bDeg=bDegI || (aJD1==j || aJD2==j);
340 aIP.SetDegenerated(bDeg);
347 TPoints.SetNbItems(iCnt);
349 IntCurveSurface_ThePolyhedronOfHInter polyhedron(aS, Upars, Vpars);
351 aTol=polyhedron.DeflectionOverEstimation();
354 Standard_Real a1,a2,a3,b1,b2,b3;
356 aBox.Get(a1,a2,a3,b1,b2,b3);
357 aBox.Update(a1-aTol,a2-aTol,a3-aTol,b1+aTol,b2+aTol,b3+aTol);
358 aBox.Enlarge(MyTolerance);
361 //=======================================================================
362 //function : FillArrayOfPnt
363 //purpose : Compute points on one surface and fill an array of points
364 // REMPLISSAGE DU TABLEAU DE POINTS
365 //=======================================================================
366 void IntPolyh_MaillageAffinage::FillArrayOfPnt
367 (const Standard_Integer SurfID,
368 const Standard_Boolean isShiftFwd,
369 const TColStd_Array1OfReal& Upars,
370 const TColStd_Array1OfReal& Vpars)
372 Standard_Boolean bDegI, bDeg;
373 Standard_Integer aNbU, aNbV, iCnt, i, j;
374 Standard_Integer aID1, aID2, aJD1, aJD2;
375 Standard_Real Tol, resol, u0, v0, u1, v1, aU, aV, aMag;
376 Standard_Real aX, aY, aZ;
378 gp_Vec aDU, aDV, aNorm;
380 aNbU=(SurfID==1)? NbSamplesU1:NbSamplesU2;
381 aNbV=(SurfID==1)? NbSamplesV1:NbSamplesV2;
382 Bnd_Box& aBox = (SurfID==1) ? MyBox1 : MyBox2;
383 Handle(Adaptor3d_HSurface) aS=(SurfID==1)? MaSurface1:MaSurface2;
384 IntPolyh_ArrayOfPoints &TPoints=(SurfID==1)? TPoints1:TPoints2;
386 resol = gp::Resolution();
392 IntCurveSurface_ThePolyhedronOfHInter polyhedron(aS, Upars, Vpars);
393 Tol=polyhedron.DeflectionOverEstimation();
398 DegeneratedIndex(Vpars, aNbV, aS, 1, aJD1, aJD2);
399 if (!(aJD1 || aJD2)) {
400 DegeneratedIndex(Upars, aNbU, aS, 2, aID1, aID2);
403 TPoints.Init(aNbU*aNbV);
405 for(i=1; i<=aNbU; ++i){
406 bDegI=(aID1==i || aID2==i);
408 for(j=1; j<=aNbV; ++j){
410 aS->D1(aU, aV, aP, aDU, aDV);
412 aNorm = aDU.Crossed(aDV);
413 aMag = aNorm.Magnitude();
416 aNorm.Multiply(Tol*1.5);
422 aP.Translate(aNorm.Reversed());
426 IntPolyh_Point& aIP=TPoints[iCnt];
427 aP.Coord(aX, aY, aZ);
428 aIP.Set(aX, aY, aZ, aU, aV);
430 bDeg=bDegI || (aJD1==j || aJD2==j);
432 aIP.SetDegenerated(bDeg);
439 TPoints.SetNbItems(iCnt);
443 Standard_Real a1,a2,a3,b1,b2,b3;
445 aBox.Get(a1,a2,a3,b1,b2,b3);
446 aBox.Update(a1-Tol,a2-Tol,a3-Tol,b1+Tol,b2+Tol,b3+Tol);
447 aBox.Enlarge(MyTolerance);
449 //=======================================================================
450 //function : CommonBox
451 //purpose : Compute the common box witch is the intersection
452 // of the two bounding boxes, and mark the points of
453 // the two surfaces that are inside.
454 // REJECTION BOUNDING BOXES
455 // DETERMINATION OF THE COMMON BOX
456 //=======================================================================
457 void IntPolyh_MaillageAffinage::CommonBox (const Bnd_Box &,
466 Standard_Real x10,y10,z10,x11,y11,z11;
467 Standard_Real x20,y20,z20,x21,y21,z21;
469 MyBox1.Get(x10,y10,z10,x11,y11,z11);
470 MyBox2.Get(x20,y20,z20,x21,y21,z21);
478 if((x10>x21)||(x20>x11)||(y10>y21)||(y20>y11)||(z10>z21)||(z20>z11)) {
481 if(x11<=x21) XMax=x11; else { if(x21<=x11) XMax=x21;}
482 if(x20<=x10) XMin=x10; else { if(x10<=x20) XMin=x20;}
483 if(y11<=y21) YMax=y11; else { if(y21<=y11) YMax=y21;}
484 if(y20<=y10) YMin=y10; else { if(y10<=y20) YMin=y20;}
485 if(z11<=z21) ZMax=z11; else { if(z21<=z11) ZMax=z21;}
486 if(z20<=z10) ZMin=z10; else { if(z10<=z20) ZMin=z20;}
488 if(((XMin==XMax)&&(!(YMin==YMax)&&!(ZMin==ZMax)))
489 ||((YMin==YMax)&&(!(XMin==XMax)&&!(ZMin==ZMax)))//ou exclusif ??
490 ||((ZMin==ZMax)&&(!(XMin==XMax)&&!(YMin==YMax)))) {
498 //extension of the box
499 if( (X==0)&&(Y!=0) ) X=Y*0.1;
500 else if( (X==0)&&(Z!=0) ) X=Z*0.1;
503 if( (Y==0)&&(X!=0) ) Y=X*0.1;
504 else if( (Y==0)&&(Z!=0) ) Y=Z*0.1;
507 if( (Z==0)&&(X!=0) ) Z=X*0.1;
508 else if( (Z==0)&&(Y!=0) ) Z=Y*0.1;
512 if( (X==0)&&(Y==0)&&(Z==0) ) {
520 //Marking of points included in the common
521 const Standard_Integer FinTP1 = TPoints1.NbItems();
522 // for(Standard_Integer i=0; i<FinTP1; i++) {
524 for( i=0; i<FinTP1; i++) {
525 IntPolyh_Point & Pt1 = TPoints1[i];
552 Pt1.SetPartOfCommon(r);
555 const Standard_Integer FinTP2 = TPoints2.NbItems();
556 for(Standard_Integer ii=0; ii<FinTP2; ii++) {
557 IntPolyh_Point & Pt2 = TPoints2[ii];
585 Pt2.SetPartOfCommon(rr);
588 //=======================================================================
589 //function : FillArrayOfEdges
590 //purpose : Compute edges from the array of points
591 // FILL THE ARRAY OF EDGES
592 //=======================================================================
593 void IntPolyh_MaillageAffinage::FillArrayOfEdges
594 (const Standard_Integer SurfID)
597 IntPolyh_ArrayOfEdges &TEdges=(SurfID==1)? TEdges1:TEdges2;
598 Standard_Integer NbSamplesU=(SurfID==1)? NbSamplesU1:NbSamplesU2;
599 Standard_Integer NbSamplesV=(SurfID==1)? NbSamplesV1:NbSamplesV2;
601 //NbEdges = 3 + 3*(NbSamplesV-2) + 3*(NbSamplesU-2) +
602 // + 3*(NbSamplesU-2)*(NbSamplesV-2) + (NbSamplesV-1) + (NbSamplesU-1);
603 //NbSamplesU and NbSamples cannot be less than 2, so
604 Standard_Integer NbEdges = 3*NbSamplesU*NbSamplesV - 2*(NbSamplesU+NbSamplesV) + 1;
605 TEdges.Init(NbEdges);
607 Standard_Integer CpteurTabEdges=0;
610 TEdges[CpteurTabEdges].SetFirstPoint(0); // U V
611 TEdges[CpteurTabEdges].SetSecondPoint(1); // U V+1
612 // TEdges[CpteurTabEdges].SetFirstTriangle(-1);
613 TEdges[CpteurTabEdges].SetSecondTriangle(0);
616 TEdges[CpteurTabEdges].SetFirstPoint(0); // U V
617 TEdges[CpteurTabEdges].SetSecondPoint(NbSamplesV); // U+1 V
618 TEdges[CpteurTabEdges].SetFirstTriangle(0);
619 TEdges[CpteurTabEdges].SetSecondTriangle(1);
622 TEdges[CpteurTabEdges].SetFirstPoint(0); // U V
623 TEdges[CpteurTabEdges].SetSecondPoint(NbSamplesV+1); // U+1 V+1
624 TEdges[CpteurTabEdges].SetFirstTriangle(1);
625 // TEdges[CpteurTabEdges].SetSecondTriangle(-1);
629 Standard_Integer PntInit=1;
630 Standard_Integer BoucleMeshV;
631 for(BoucleMeshV=1; BoucleMeshV<NbSamplesV-1;BoucleMeshV++){
632 TEdges[CpteurTabEdges].SetFirstPoint(PntInit); // U V
633 TEdges[CpteurTabEdges].SetSecondPoint(PntInit+1); // U V+1
634 // TEdges[CpteurTabEdges].SetFirstTriangle(-1);
635 TEdges[CpteurTabEdges].SetSecondTriangle(BoucleMeshV*2);
638 TEdges[CpteurTabEdges].SetFirstPoint(PntInit); // U V
639 TEdges[CpteurTabEdges].SetSecondPoint(PntInit+NbSamplesV+1); // U+1 V+1
640 TEdges[CpteurTabEdges].SetFirstTriangle(BoucleMeshV*2);
641 TEdges[CpteurTabEdges].SetSecondTriangle(BoucleMeshV*2+1);
644 TEdges[CpteurTabEdges].SetFirstPoint(PntInit); // U V
645 TEdges[CpteurTabEdges].SetSecondPoint(PntInit+NbSamplesV); // U+1 V
646 TEdges[CpteurTabEdges].SetFirstTriangle(BoucleMeshV*2+1);
647 TEdges[CpteurTabEdges].SetSecondTriangle(BoucleMeshV*2-2);
654 for(BoucleMeshV=1; BoucleMeshV<NbSamplesU-1;BoucleMeshV++){
655 TEdges[CpteurTabEdges].SetFirstPoint(PntInit); // U V
656 TEdges[CpteurTabEdges].SetSecondPoint(PntInit+1); // U V+1
657 TEdges[CpteurTabEdges].SetFirstTriangle((BoucleMeshV-1)*(NbSamplesV-1)*2+1);
658 TEdges[CpteurTabEdges].SetSecondTriangle(BoucleMeshV*(NbSamplesV-1)*2);
661 TEdges[CpteurTabEdges].SetFirstPoint(PntInit); // U V
662 TEdges[CpteurTabEdges].SetSecondPoint(PntInit+NbSamplesV+1); // U+1 V+1
663 TEdges[CpteurTabEdges].SetFirstTriangle(BoucleMeshV*(NbSamplesV-1)*2);
664 TEdges[CpteurTabEdges].SetSecondTriangle(BoucleMeshV*(NbSamplesV-1)*2+1);
667 TEdges[CpteurTabEdges].SetFirstPoint(PntInit); // U V
668 TEdges[CpteurTabEdges].SetSecondPoint(PntInit+NbSamplesV); // U+1 V
669 TEdges[CpteurTabEdges].SetFirstTriangle(BoucleMeshV*(NbSamplesV-1)*2+1);
670 // TEdges[CpteurTabEdges].SetSecondTriangle(-1);
675 PntInit=NbSamplesV+1;
676 //To provide recursion I associate a point with three edges
677 for(Standard_Integer BoucleMeshU=1; BoucleMeshU<NbSamplesU-1; BoucleMeshU++){
678 for(Standard_Integer BoucleMeshV=1; BoucleMeshV<NbSamplesV-1;BoucleMeshV++){
679 TEdges[CpteurTabEdges].SetFirstPoint(PntInit); // U V
680 TEdges[CpteurTabEdges].SetSecondPoint(PntInit+1); // U V+1
681 TEdges[CpteurTabEdges].SetFirstTriangle((NbSamplesV-1)*2*(BoucleMeshU-1)+BoucleMeshV*2+1);
682 TEdges[CpteurTabEdges].SetSecondTriangle((NbSamplesV-1)*2*BoucleMeshU+BoucleMeshV*2);
685 TEdges[CpteurTabEdges].SetFirstPoint(PntInit); // U V
686 TEdges[CpteurTabEdges].SetSecondPoint(PntInit+NbSamplesV+1); // U+1 V+1
687 TEdges[CpteurTabEdges].SetFirstTriangle((NbSamplesV-1)*2*BoucleMeshU+BoucleMeshV*2);
688 TEdges[CpteurTabEdges].SetSecondTriangle((NbSamplesV-1)*2*BoucleMeshU+BoucleMeshV*2+1);
691 TEdges[CpteurTabEdges].SetFirstPoint(PntInit); // U V
692 TEdges[CpteurTabEdges].SetSecondPoint(PntInit+NbSamplesV); // U+1 V
693 TEdges[CpteurTabEdges].SetFirstTriangle((NbSamplesV-1)*2*BoucleMeshU+BoucleMeshV*2+1);
694 TEdges[CpteurTabEdges].SetSecondTriangle((NbSamplesV-1)*2*BoucleMeshU+BoucleMeshV*2-2);
696 PntInit++;//Pass to the next point
698 PntInit++;//Pass the last point of the column
699 PntInit++;//Pass the first point of the next column
703 PntInit=(NbSamplesU-1)*NbSamplesV; //point U=u1 V=0
704 for(BoucleMeshV=0; BoucleMeshV<NbSamplesV-1; BoucleMeshV++){
705 TEdges[CpteurTabEdges].SetFirstPoint(PntInit); //U=u1 V
706 TEdges[CpteurTabEdges].SetSecondPoint(PntInit+1); //U=u1 V+1
707 TEdges[CpteurTabEdges].SetFirstTriangle((NbSamplesU-2)*(NbSamplesV-1)*2+BoucleMeshV*2+1);
708 // TEdges[CpteurTabEdges].SetSecondTriangle(-1);
714 for(BoucleMeshV=0; BoucleMeshV<NbSamplesU-1;BoucleMeshV++){
715 TEdges[CpteurTabEdges].SetFirstPoint(NbSamplesV-1+BoucleMeshV*NbSamplesV); // U V=v1
716 TEdges[CpteurTabEdges].SetSecondPoint(NbSamplesV-1+(BoucleMeshV+1)*NbSamplesV); //U+1 V=v1
717 // TEdges[CpteurTabEdges].SetFirstTriangle(-1);
718 TEdges[CpteurTabEdges].SetSecondTriangle(BoucleMeshV*2*(NbSamplesV-1)+(NbSamplesV-2)*2);
721 TEdges.SetNbItems(CpteurTabEdges);
724 //=======================================================================
725 //function : FillArrayOfTriangles
726 //purpose : Compute triangles from the array of points, and --
727 // mark the triangles that use marked points by the
728 // CommonBox function.
729 // FILL THE ARRAY OF TRIANGLES
730 //=======================================================================
731 void IntPolyh_MaillageAffinage::FillArrayOfTriangles
732 (const Standard_Integer SurfID)
734 Standard_Integer CpteurTabTriangles=0;
735 Standard_Integer PntInit=0;
737 IntPolyh_ArrayOfPoints &TPoints=(SurfID==1)? TPoints1:TPoints2;
738 IntPolyh_ArrayOfTriangles &TTriangles=(SurfID==1)? TTriangles1:TTriangles2;
739 Standard_Integer NbSamplesU=(SurfID==1)? NbSamplesU1:NbSamplesU2;
740 Standard_Integer NbSamplesV=(SurfID==1)? NbSamplesV1:NbSamplesV2;
742 TTriangles.Init(2*(NbSamplesU-1)*(NbSamplesV-1));
743 //To provide recursion, I associate a point with two triangles
744 for(Standard_Integer BoucleMeshU=0; BoucleMeshU<NbSamplesU-1; BoucleMeshU++){
745 for(Standard_Integer BoucleMeshV=0; BoucleMeshV<NbSamplesV-1;BoucleMeshV++){
748 TTriangles[CpteurTabTriangles].SetFirstPoint(PntInit); // U V
749 TTriangles[CpteurTabTriangles].SetSecondPoint(PntInit+1); // U V+1
750 TTriangles[CpteurTabTriangles].SetThirdPoint(PntInit+NbSamplesV+1); // U+1 V+1
752 // IF ITS EDGE CONTACTS WITH THE COMMON BOX IP REMAINS = A 1
753 if( ( (TPoints[PntInit].PartOfCommon()) & (TPoints[PntInit+1].PartOfCommon()) )
754 &&( (TPoints[PntInit+1].PartOfCommon()) & (TPoints[PntInit+NbSamplesV+1].PartOfCommon()))
755 &&( (TPoints[PntInit+NbSamplesV+1].PartOfCommon()) & (TPoints[PntInit].PartOfCommon())) )
757 TTriangles[CpteurTabTriangles].SetIndiceIntersectionPossible(0);
759 CpteurTabTriangles++;
762 TTriangles[CpteurTabTriangles].SetFirstPoint(PntInit); // U V
763 TTriangles[CpteurTabTriangles].SetSecondPoint(PntInit+NbSamplesV+1); // U+1 V+1
764 TTriangles[CpteurTabTriangles].SetThirdPoint(PntInit+NbSamplesV); // U+1 V
767 if( ( (TPoints[PntInit].PartOfCommon()) & (TPoints[PntInit+NbSamplesV+1].PartOfCommon()) )
768 &&( (TPoints[PntInit+NbSamplesV+1].PartOfCommon()) & (TPoints[PntInit+NbSamplesV].PartOfCommon()))
769 &&( (TPoints[PntInit+NbSamplesV].PartOfCommon()) & (TPoints[PntInit].PartOfCommon())) )
770 TTriangles[CpteurTabTriangles].SetIndiceIntersectionPossible(0);
773 CpteurTabTriangles++;
775 PntInit++;//Pass to the next point
777 PntInit++;//Pass the last point of the column
779 TTriangles.SetNbItems(CpteurTabTriangles);
780 const Standard_Integer FinTT = TTriangles.NbItems();
784 //=======================================================================
785 //function : LinkEdges2Triangles
786 //purpose : fill the edge fields in Triangle object for the
787 // two array of triangles.
788 //=======================================================================
789 void IntPolyh_MaillageAffinage::LinkEdges2Triangles()
791 const Standard_Integer FinTT1 = TTriangles1.NbItems();
792 const Standard_Integer FinTT2 = TTriangles2.NbItems();
794 for(Standard_Integer uiui1=0; uiui1<FinTT1; uiui1++) {
795 IntPolyh_Triangle & MyTriangle1=TTriangles1[uiui1];
796 if ( (MyTriangle1.FirstEdge()) == -1 ) {
797 MyTriangle1.SetEdgeandOrientation(1,TEdges1);
798 MyTriangle1.SetEdgeandOrientation(2,TEdges1);
799 MyTriangle1.SetEdgeandOrientation(3,TEdges1);
802 for(Standard_Integer uiui2=0; uiui2<FinTT2; uiui2++) {
803 IntPolyh_Triangle & MyTriangle2=TTriangles2[uiui2];
804 if ( (MyTriangle2.FirstEdge()) == -1 ) {
805 MyTriangle2.SetEdgeandOrientation(1,TEdges2);
806 MyTriangle2.SetEdgeandOrientation(2,TEdges2);
807 MyTriangle2.SetEdgeandOrientation(3,TEdges2);
811 //=======================================================================
812 //function : CommonPartRefinement
813 //purpose : Refine systematicaly all marked triangles of both surfaces
814 // REFINING OF THE COMMON
815 //=======================================================================
816 void IntPolyh_MaillageAffinage::CommonPartRefinement()
818 Standard_Integer FinInit1 = TTriangles1.NbItems();
819 for(Standard_Integer i=0; i<FinInit1; i++) {
820 if(TTriangles1[i].IndiceIntersectionPossible()!=0)
821 TTriangles1[i].MiddleRefinement(i,MaSurface1,TPoints1,TTriangles1,TEdges1);
824 Standard_Integer FinInit2=TTriangles2.NbItems();
825 for(Standard_Integer ii=0; ii<FinInit2; ii++) {
826 if(TTriangles2[ii].IndiceIntersectionPossible()!=0)
827 TTriangles2[ii].MiddleRefinement(ii,MaSurface2,TPoints2,TTriangles2,TEdges2);
831 //=======================================================================
832 //function : LocalSurfaceRefinement
833 //purpose : Refine systematicaly all marked triangles of ONE surface
834 //=======================================================================
835 void IntPolyh_MaillageAffinage::LocalSurfaceRefinement(const Standard_Integer SurfID) {
836 //refine locally, but systematically the chosen surface
838 const Standard_Integer FinInit1 = TTriangles1.NbItems();
839 for(Standard_Integer i=0; i<FinInit1; i++) {
840 if(TTriangles1[i].IndiceIntersectionPossible()!=0)
841 TTriangles1[i].MiddleRefinement(i,MaSurface1,TPoints1,TTriangles1,TEdges1);
846 const Standard_Integer FinInit2 = TTriangles2.NbItems();
847 for(Standard_Integer ii=0; ii<FinInit2; ii++) {
848 if(TTriangles2[ii].IndiceIntersectionPossible()!=0)
849 TTriangles2[ii].MiddleRefinement(ii,MaSurface2,TPoints2,TTriangles2,TEdges2);
853 //=======================================================================
854 //function : ComputeDeflections
855 //purpose : Compute deflection for all triangles of one
856 // surface,and sort min and max of deflections
858 // Calculation of the deflection of all triangles
859 // --> deflection max
860 // --> deflection min
861 //=======================================================================
862 void IntPolyh_MaillageAffinage::ComputeDeflections
863 (const Standard_Integer SurfID)
865 Handle(Adaptor3d_HSurface) MaSurface=(SurfID==1)? MaSurface1:MaSurface2;
866 IntPolyh_ArrayOfPoints &TPoints=(SurfID==1)? TPoints1:TPoints2;
867 IntPolyh_ArrayOfTriangles &TTriangles=(SurfID==1)? TTriangles1:TTriangles2;
868 Standard_Real &FlecheMin=(SurfID==1)? FlecheMin1:FlecheMin2;
869 Standard_Real &FlecheMoy=(SurfID==1)? FlecheMoy1:FlecheMoy2;
870 Standard_Real &FlecheMax=(SurfID==1)? FlecheMax1:FlecheMax2;
872 Standard_Integer CpteurTabFleche=0;
873 FlecheMax=-RealLast();
874 FlecheMin=RealLast();
876 const Standard_Integer FinTT = TTriangles.NbItems();
878 for(CpteurTabFleche=0; CpteurTabFleche<FinTT; CpteurTabFleche++) {
879 IntPolyh_Triangle &Triangle = TTriangles[CpteurTabFleche];
880 if ( Triangle.GetFleche() < 0) { //pas normal
884 Triangle.TriangleDeflection(MaSurface, TPoints);
885 Standard_Real Fleche=Triangle.GetFleche();
887 if (Fleche > FlecheMax)
889 if (Fleche < FlecheMin)
894 //=======================================================================
895 //function : TrianglesDeflectionsRefinementBSB
896 //purpose : Refine both surfaces using BoundSortBox as --
897 // rejection. The criterions used to refine a --
898 // triangle are: The deflection The size of the --
899 // bounding boxes (one surface may be very small
900 // compared to the other)
901 //=======================================================================
902 void IntPolyh_MaillageAffinage::TrianglesDeflectionsRefinementBSB()
904 const Standard_Integer FinTT1 = TTriangles1.NbItems();
905 const Standard_Integer FinTT2 = TTriangles2.NbItems();
907 ComputeDeflections(1);
908 // To estimate a surface in general it can be interesting
909 //to calculate all deflections
910 //-- Check deflection at output
912 Standard_Real FlecheCritique1;
913 if(FlecheMin1>FlecheMax1) {
916 else {//fleche min + (flechemax-flechemin) * 80/100
917 FlecheCritique1 = FlecheMin1*0.2+FlecheMax1*0.8;
920 ComputeDeflections(2);
921 //-- Check arrows at output
923 Standard_Real FlecheCritique2;
924 if(FlecheMin2>FlecheMax2) {
928 else {//fleche min + (flechemax-flechemin) * 80/100
929 FlecheCritique2 = FlecheMin2*0.2+FlecheMax2*0.8;
933 Bnd_BoundSortBox BndBSB;
934 Standard_Real diag1,diag2;
935 Standard_Real x0,y0,z0,x1,y1,z1;
937 //The greatest of two bounding boxes created in FillArrayOfPoints is found.
938 //Then this value is weighted depending on the discretization
939 //(NbSamplesU and NbSamplesV)
940 MyBox1.Get(x0,y0,z0,x1,y1,z1);
941 x0-=x1; y0-=y1; z0-=z1;
942 diag1=x0*x0+y0*y0+z0*z0;
943 const Standard_Real NbSamplesUV1=Standard_Real(NbSamplesU1) * Standard_Real(NbSamplesV1);
946 MyBox2.Get(x0,y0,z0,x1,y1,z1);
947 x0-=x1; y0-=y1; z0-=z1;
948 diag2=x0*x0+y0*y0+z0*z0;
949 const Standard_Real NbSamplesUV2=Standard_Real(NbSamplesU2) * Standard_Real(NbSamplesV2);
952 //-- The surface with the greatest bounding box is "discretized"
954 //Standard_Integer NbInterTentees=0;
958 if(FlecheCritique2<diag1) {//the corresponding sizes are not too disproportional
960 Handle(Bnd_HArray1OfBox) HBnd = new Bnd_HArray1OfBox(1,FinTT2);
962 for(Standard_Integer i=0; i<FinTT2; i++){
963 if (TTriangles2[i].IndiceIntersectionPossible()!=0) {
965 const IntPolyh_Triangle& T=TTriangles2[i];
966 const IntPolyh_Point& PA=TPoints2[T.FirstPoint()];
967 const IntPolyh_Point& PB=TPoints2[T.SecondPoint()];
968 const IntPolyh_Point& PC=TPoints2[T.ThirdPoint()];
969 gp_Pnt pntA(PA.X(),PA.Y(),PA.Z());
970 gp_Pnt pntB(PB.X(),PB.Y(),PB.Z());
971 gp_Pnt pntC(PC.X(),PC.Y(),PC.Z());
972 b.Add(pntA);//Box b, which contains triangle i of surface 2 is created./
975 b.Enlarge(T.GetFleche()+MyTolerance);
976 HBnd->SetValue(i+1,b);//Box b is added in the array HBnd
980 //Inititalization of the boundary, sorting of boxes
981 BndBSB.Initialize(HBnd);//contains boxes of 2
983 Standard_Integer FinTT1Init=FinTT1;
984 for(Standard_Integer i_S1=0; i_S1<FinTT1Init; i_S1++) {
985 if(TTriangles1[i_S1].IndiceIntersectionPossible()!=0) {
986 //-- Loop on the boxes of mesh 1
988 const IntPolyh_Triangle& T=TTriangles1[i_S1];
989 const IntPolyh_Point& PA=TPoints1[T.FirstPoint()];
990 const IntPolyh_Point& PB=TPoints1[T.SecondPoint()];
991 const IntPolyh_Point& PC=TPoints1[T.ThirdPoint()];
992 gp_Pnt pntA(PA.X(),PA.Y(),PA.Z());
993 gp_Pnt pntB(PB.X(),PB.Y(),PB.Z());
994 gp_Pnt pntC(PC.X(),PC.Y(),PC.Z());
998 b.Enlarge(T.GetFleche());
999 //-- List of boxes of 2, which touch this box (of 1)
1000 const TColStd_ListOfInteger& ListeOf2 = BndBSB.Compare(b);
1002 if((ListeOf2.IsEmpty())==0) {
1003 IntPolyh_Triangle &Triangle1 = TTriangles1[i_S1];
1004 if(Triangle1.GetFleche()>FlecheCritique1)
1005 Triangle1.MiddleRefinement(i_S1,MaSurface1,TPoints1,
1006 TTriangles1, TEdges1);
1008 for (TColStd_ListIteratorOfListOfInteger Iter(ListeOf2);
1011 Standard_Integer i_S2=Iter.Value()-1;
1012 //if the box of s1 contacts with the boxes of s2
1013 //the arrow of the triangle is checked
1014 IntPolyh_Triangle & Triangle2 = TTriangles2[i_S2];
1015 if(Triangle2.IndiceIntersectionPossible()!=0)
1016 if(Triangle2.GetFleche()>FlecheCritique2)
1017 Triangle2.MiddleRefinement( i_S2, MaSurface2, TPoints2,
1018 TTriangles2, TEdges2);
1025 //--------------------------------------------------------------------
1026 //FlecheCritique2 > diag1
1030 Handle(Bnd_HArray1OfBox) HBnd = new Bnd_HArray1OfBox(1,FinTT2);
1032 for(Standard_Integer i=0; i<FinTT2; i++){
1033 if (TTriangles2[i].IndiceIntersectionPossible()!=0) {
1035 const IntPolyh_Triangle& T=TTriangles2[i];
1036 const IntPolyh_Point& PA=TPoints2[T.FirstPoint()];
1037 const IntPolyh_Point& PB=TPoints2[T.SecondPoint()];
1038 const IntPolyh_Point& PC=TPoints2[T.ThirdPoint()];
1039 gp_Pnt pntA(PA.X(),PA.Y(),PA.Z());
1040 gp_Pnt pntB(PB.X(),PB.Y(),PB.Z());
1041 gp_Pnt pntC(PC.X(),PC.Y(),PC.Z());
1042 b.Add(pntA);//Box b, which contains triangle i of surface 2 is created/
1045 b.Enlarge(T.GetFleche()+MyTolerance);
1046 //-- BndBSB.Add(b,i+1);
1047 HBnd->SetValue(i+1,b);//Box b is added in array HBnd
1051 //Inititalization of the ouput bounding box
1052 BndBSB.Initialize(HBnd);//contains boxes of 2
1055 //The bounding box Be1 of surface1 is compared BSB of surface2
1056 const TColStd_ListOfInteger& ListeOf2 = BndBSB.Compare(MyBox1);
1058 if((ListeOf2.IsEmpty())==0) {
1059 //if the bounding box Be1 of s1 contacts with
1060 //the boxes of s2 the deflection of triangle of s2 is checked
1062 // Be1 is very small in relation to Be2
1063 //The criterion of refining for surface2 depends on the size of Be1
1064 //As it is known that this criterion should be minimized,
1065 //the smallest side of the bounding box is taken
1066 Standard_Real x0,x1,y0,y1,z0,z1;
1067 MyBox1.Get(x0,y0,z0,x1,y1,z1);
1068 Standard_Real dx=Abs(x1-x0);
1069 Standard_Real dy=Abs(y1-y0);
1070 Standard_Real diag=Abs(z1-z0);
1071 Standard_Real dd=-1.0;
1076 if (diag>dd) diag=dd;
1078 //if Be1 contacts with the boxes of s2, the deflection
1079 //of the triangles of s2 is checked (greater)
1080 //in relation to the size of Be1 (smaller)
1081 for (TColStd_ListIteratorOfListOfInteger Iter(ListeOf2);
1084 Standard_Integer i_S2=Iter.Value()-1;
1086 IntPolyh_Triangle & Triangle2=TTriangles2[i_S2];
1087 if(Triangle2.IndiceIntersectionPossible()) {
1089 //calculation of the criterion of refining
1090 //The deflection of the greater is compared to the size of the smaller
1091 Standard_Real CritereAffinage=0.0;
1092 Standard_Real DiagPonderation=0.5;
1093 CritereAffinage = diag*DiagPonderation;
1094 if(Triangle2.GetFleche()>CritereAffinage)
1095 Triangle2.MultipleMiddleRefinement2(CritereAffinage, MyBox1, i_S2,
1096 MaSurface2, TPoints2,
1097 TTriangles2,TEdges2);
1099 else Triangle2.MiddleRefinement(i_S2,MaSurface2,TPoints2,
1100 TTriangles2, TEdges2);
1108 else { //-- The greater is discretised
1110 if(FlecheCritique1<diag2) {//the respective sizes are not to much disproportional
1112 Handle(Bnd_HArray1OfBox) HBnd = new Bnd_HArray1OfBox(1,FinTT1);
1114 for(Standard_Integer i=0; i<FinTT1; i++){
1115 if(TTriangles1[i].IndiceIntersectionPossible()!=0) {
1117 const IntPolyh_Triangle& T=TTriangles1[i];
1118 const IntPolyh_Point& PA=TPoints1[T.FirstPoint()];
1119 const IntPolyh_Point& PB=TPoints1[T.SecondPoint()];
1120 const IntPolyh_Point& PC=TPoints1[T.ThirdPoint()];
1121 gp_Pnt pntA(PA.X(),PA.Y(),PA.Z());
1122 gp_Pnt pntB(PB.X(),PB.Y(),PB.Z());
1123 gp_Pnt pntC(PC.X(),PC.Y(),PC.Z());
1124 b.Add(pntA);//Box b, which contains triangle i of surface 2 is created.
1127 b.Enlarge(T.GetFleche()+MyTolerance);
1128 HBnd->SetValue(i+1,b);//Boite b is added in the array HBnd
1131 BndBSB.Initialize(HBnd);
1133 Standard_Integer FinTT2init=FinTT2;
1134 for(Standard_Integer i_S2=0; i_S2<FinTT2init; i_S2++) {
1135 if (TTriangles2[i_S2].IndiceIntersectionPossible()!=0) {
1136 //-- Loop on the boxes of mesh 2
1138 const IntPolyh_Triangle& T=TTriangles2[i_S2];
1139 const IntPolyh_Point& PA=TPoints2[T.FirstPoint()];
1140 const IntPolyh_Point& PB=TPoints2[T.SecondPoint()];
1141 const IntPolyh_Point& PC=TPoints2[T.ThirdPoint()];
1142 gp_Pnt pntA(PA.X(),PA.Y(),PA.Z());
1143 gp_Pnt pntB(PB.X(),PB.Y(),PB.Z());
1144 gp_Pnt pntC(PC.X(),PC.Y(),PC.Z());
1148 b.Enlarge(T.GetFleche()+MyTolerance);
1149 //-- List of boxes of 1 touching this box (of 2)
1150 const TColStd_ListOfInteger& ListeOf1 = BndBSB.Compare(b);
1151 IntPolyh_Triangle & Triangle2=TTriangles2[i_S2];
1152 if((ListeOf1.IsEmpty())==0) {
1154 if(Triangle2.GetFleche()>FlecheCritique2)
1155 Triangle2.MiddleRefinement(i_S2,MaSurface2,TPoints2,
1156 TTriangles2, TEdges2);
1158 for (TColStd_ListIteratorOfListOfInteger Iter(ListeOf1);
1161 Standard_Integer i_S1=Iter.Value()-1;
1162 IntPolyh_Triangle & Triangle1=TTriangles1[i_S1];
1163 if (Triangle1.IndiceIntersectionPossible())
1164 if(Triangle1.GetFleche()>FlecheCritique1)
1165 Triangle1.MiddleRefinement(i_S1,MaSurface1,TPoints1,
1166 TTriangles1, TEdges1);
1172 //-----------------------------------------------------------------------------
1173 else {// FlecheCritique1>diag2
1176 Handle(Bnd_HArray1OfBox) HBnd = new Bnd_HArray1OfBox(1,FinTT1);
1178 for(Standard_Integer i=0; i<FinTT1; i++){
1179 if (TTriangles1[i].IndiceIntersectionPossible()!=0) {
1181 const IntPolyh_Triangle& T=TTriangles1[i];
1182 const IntPolyh_Point& PA=TPoints1[T.FirstPoint()];
1183 const IntPolyh_Point& PB=TPoints1[T.SecondPoint()];
1184 const IntPolyh_Point& PC=TPoints1[T.ThirdPoint()];
1185 gp_Pnt pntA(PA.X(),PA.Y(),PA.Z());
1186 gp_Pnt pntB(PB.X(),PB.Y(),PB.Z());
1187 gp_Pnt pntC(PC.X(),PC.Y(),PC.Z());
1188 b.Add(pntA);//Box b, which contains triangle i of surface 1 is created./
1191 b.Enlarge(T.GetFleche()+MyTolerance);
1192 HBnd->SetValue(i+1,b);//Box b is added in the array HBnd
1196 //Inititalisation of the boundary output box
1197 BndBSB.Initialize(HBnd);//contains boxes of 1
1199 //Bounding box Be2 of surface2 is compared to BSB of surface1
1200 const TColStd_ListOfInteger& ListeOf1 = BndBSB.Compare(MyBox2);
1202 if((ListeOf1.IsEmpty())==0) {
1203 //if the bounding box Be2 of s2 contacts
1204 //with boxes of s1 the deflection of the triangle of s1 is checked
1206 // Be2 is very small compared to Be1
1207 //The criterion of refining for surface1 depends on the size of Be2
1208 //As this criterion should be minimized,
1209 //the smallest side of the bounding box is taken
1210 Standard_Real x0,x1,y0,y1,z0,z1;
1211 MyBox2.Get(x0,y0,z0,x1,y1,z1);
1212 Standard_Real dx=Abs(x1-x0);
1213 Standard_Real dy=Abs(y1-y0);
1214 Standard_Real diag=Abs(z1-z0);
1215 Standard_Real dd=-1.0;
1220 if (diag>dd) diag=dd;
1222 //if Be2 contacts with boxes of s1, the deflection of
1223 //triangles of s1 (greater) is checked
1224 //comparatively to the size of Be2 (smaller).
1225 for (TColStd_ListIteratorOfListOfInteger Iter(ListeOf1);
1228 Standard_Integer i_S1=Iter.Value()-1;
1230 IntPolyh_Triangle & Triangle1=TTriangles1[i_S1];
1231 if(Triangle1.IndiceIntersectionPossible()) {
1233 //calculation of the criterion of refining
1234 //The deflection of the greater is compared
1235 //with the size of the smaller.
1236 Standard_Real CritereAffinage=0.0;
1237 Standard_Real DiagPonderation=0.5;
1238 CritereAffinage = diag*DiagPonderation;;
1239 if(Triangle1.GetFleche()>CritereAffinage)
1240 Triangle1.MultipleMiddleRefinement2(CritereAffinage,MyBox2, i_S1,
1241 MaSurface1, TPoints1,
1242 TTriangles1, TEdges1);
1244 else Triangle1.MiddleRefinement(i_S1,MaSurface1,TPoints1,
1245 TTriangles1, TEdges1);
1253 //=======================================================================
1255 //purpose : This function is used for the function project6
1256 //=======================================================================
1257 inline Standard_Real maxSR(const Standard_Real a,
1258 const Standard_Real b,
1259 const Standard_Real c)
1261 Standard_Real t = a;
1266 //=======================================================================
1268 //purpose : This function is used for the function project6
1269 //=======================================================================
1270 inline Standard_Real minSR(const Standard_Real a,
1271 const Standard_Real b,
1272 const Standard_Real c)
1274 Standard_Real t = a;
1279 //=======================================================================
1280 //function : project6
1281 //purpose : This function is used for the function TriContact
1282 //=======================================================================
1283 Standard_Integer project6(const IntPolyh_Point &ax,
1284 const IntPolyh_Point &p1,
1285 const IntPolyh_Point &p2,
1286 const IntPolyh_Point &p3,
1287 const IntPolyh_Point &q1,
1288 const IntPolyh_Point &q2,
1289 const IntPolyh_Point &q3)
1291 Standard_Real P1 = ax.Dot(p1);
1292 Standard_Real P2 = ax.Dot(p2);
1293 Standard_Real P3 = ax.Dot(p3);
1294 Standard_Real Q1 = ax.Dot(q1);
1295 Standard_Real Q2 = ax.Dot(q2);
1296 Standard_Real Q3 = ax.Dot(q3);
1298 Standard_Real mx1 = maxSR(P1, P2, P3);
1299 Standard_Real mn1 = minSR(P1, P2, P3);
1300 Standard_Real mx2 = maxSR(Q1, Q2, Q3);
1301 Standard_Real mn2 = minSR(Q1, Q2, Q3);
1303 if (mn1 > mx2) return 0;
1304 if (mn2 > mx1) return 0;
1307 //=======================================================================
1308 //function : TriContact
1309 //purpose : This fonction Check if two triangles are in
1310 // contact or no, return 1 if yes, return 0
1312 //=======================================================================
1313 Standard_Integer IntPolyh_MaillageAffinage::TriContact
1314 (const IntPolyh_Point &P1,
1315 const IntPolyh_Point &P2,
1316 const IntPolyh_Point &P3,
1317 const IntPolyh_Point &Q1,
1318 const IntPolyh_Point &Q2,
1319 const IntPolyh_Point &Q3,
1320 Standard_Real &Angle) const
1323 The first triangle is (p1,p2,p3). The other is (q1,q2,q3).
1324 The edges are (e1,e2,e3) and (f1,f2,f3).
1325 The normals are n1 and m1
1326 Outwards are (g1,g2,g3) and (h1,h2,h3).*/
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);
1343 if(maxSR(P1.X(),P2.X(),P3.X())<minSR(Q1.X(),Q2.X(),Q3.X())) return(0);
1344 if(maxSR(P1.Y(),P2.Y(),P3.Y())<minSR(Q1.Y(),Q2.Y(),Q3.Y())) return(0);
1345 if(maxSR(P1.Z(),P2.Z(),P3.Z())<minSR(Q1.Z(),Q2.Z(),Q3.Z())) return(0);
1347 if(minSR(P1.X(),P2.X(),P3.X())>maxSR(Q1.X(),Q2.X(),Q3.X())) return(0);
1348 if(minSR(P1.Y(),P2.Y(),P3.Y())>maxSR(Q1.Y(),Q2.Y(),Q3.Y())) return(0);
1349 if(minSR(P1.Z(),P2.Z(),P3.Z())>maxSR(Q1.Z(),Q2.Z(),Q3.Z())) return(0);
1351 p1.SetX(P1.X() - P1.X()); p1.SetY(P1.Y() - P1.Y()); p1.SetZ(P1.Z() - P1.Z());
1352 p2.SetX(P2.X() - P1.X()); p2.SetY(P2.Y() - P1.Y()); p2.SetZ(P2.Z() - P1.Z());
1353 p3.SetX(P3.X() - P1.X()); p3.SetY(P3.Y() - P1.Y()); p3.SetZ(P3.Z() - P1.Z());
1355 q1.SetX(Q1.X() - P1.X()); q1.SetY(Q1.Y() - P1.Y()); q1.SetZ(Q1.Z() - P1.Z());
1356 q2.SetX(Q2.X() - P1.X()); q2.SetY(Q2.Y() - P1.Y()); q2.SetZ(Q2.Z() - P1.Z());
1357 q3.SetX(Q3.X() - P1.X()); q3.SetY(Q3.Y() - P1.Y()); q3.SetZ(Q3.Z() - P1.Z());
1359 e1.SetX(p2.X() - p1.X()); e1.SetY(p2.Y() - p1.Y()); e1.SetZ(p2.Z() - p1.Z());
1360 e2.SetX(p3.X() - p2.X()); e2.SetY(p3.Y() - p2.Y()); e2.SetZ(p3.Z() - p2.Z());
1361 e3.SetX(p1.X() - p3.X()); e3.SetY(p1.Y() - p3.Y()); e3.SetZ(p1.Z() - p3.Z());
1363 f1.SetX(q2.X() - q1.X()); f1.SetY(q2.Y() - q1.Y()); f1.SetZ(q2.Z() - q1.Z());
1364 f2.SetX(q3.X() - q2.X()); f2.SetY(q3.Y() - q2.Y()); f2.SetZ(q3.Z() - q2.Z());
1365 f3.SetX(q1.X() - q3.X()); f3.SetY(q1.Y() - q3.Y()); f3.SetZ(q1.Z() - q3.Z());
1367 n1.Cross(e1, e2); //normal to the first triangle
1368 m1.Cross(f1, f2); //normal to the second triangle
1387 // Now the testing is done
1389 if (!project6(n1, p1, p2, p3, q1, q2, q3)) return 0; //T2 is not higher or lower than T1
1390 if (!project6(m1, p1, p2, p3, q1, q2, q3)) return 0; //T1 is not higher of lower than T2
1392 if (!project6(ef11, p1, p2, p3, q1, q2, q3)) return 0;
1393 if (!project6(ef12, p1, p2, p3, q1, q2, q3)) return 0;
1394 if (!project6(ef13, p1, p2, p3, q1, q2, q3)) return 0;
1395 if (!project6(ef21, p1, p2, p3, q1, q2, q3)) return 0;
1396 if (!project6(ef22, p1, p2, p3, q1, q2, q3)) return 0;
1397 if (!project6(ef23, p1, p2, p3, q1, q2, q3)) return 0;
1398 if (!project6(ef31, p1, p2, p3, q1, q2, q3)) return 0;
1399 if (!project6(ef32, p1, p2, p3, q1, q2, q3)) return 0;
1400 if (!project6(ef33, p1, p2, p3, q1, q2, q3)) return 0;
1402 if (!project6(g1, p1, p2, p3, q1, q2, q3)) return 0; //T2 is outside of T1 in the plane of T1
1403 if (!project6(g2, p1, p2, p3, q1, q2, q3)) return 0; //T2 is outside of T1 in the plane of T1
1404 if (!project6(g3, p1, p2, p3, q1, q2, q3)) return 0; //T2 is outside of T1 in the plane of T1
1405 if (!project6(h1, p1, p2, p3, q1, q2, q3)) return 0; //T1 is outside of T2 in the plane of T2
1406 if (!project6(h2, p1, p2, p3, q1, q2, q3)) return 0; //T1 is outside of T2 in the plane of T2
1407 if (!project6(h3, p1, p2, p3, q1, q2, q3)) return 0; //T1 is outside of T2 in the plane of T2
1409 //Calculation of cosinus angle between two normals
1410 Standard_Real SqModn1=-1.0;
1411 Standard_Real SqModm1=-1.0;
1412 SqModn1=n1.SquareModulus();
1413 if (SqModn1>SquareMyConfusionPrecision){
1414 SqModm1=m1.SquareModulus();
1416 if (SqModm1>SquareMyConfusionPrecision) {
1417 Angle=(n1.Dot(m1))/(sqrt(SqModn1)*sqrt(SqModm1));
1421 //=======================================================================
1422 //function : TestNbPoints
1423 //purpose : This function is used by StartingPointsResearch() to control
1424 // the number of points found keep the result in conformity (1 or 2 points)
1425 // void TestNbPoints(const Standard_Integer TriSurfID,
1426 //=======================================================================
1427 void TestNbPoints(const Standard_Integer ,
1428 Standard_Integer &NbPoints,
1429 Standard_Integer &NbPointsTotal,
1430 const IntPolyh_StartPoint &Pt1,
1431 const IntPolyh_StartPoint &Pt2,
1432 IntPolyh_StartPoint &SP1,
1433 IntPolyh_StartPoint &SP2)
1435 // already checked in TriangleEdgeContact2
1436 // if( (NbPoints==2)&&(Pt1.CheckSameSP(Pt2)) ) NbPoints=1;
1442 if ( (NbPoints==1)&&(NbPointsTotal==0) ) {
1446 else if ( (NbPoints==1)&&(NbPointsTotal==1) ) {
1447 if(Pt1.CheckSameSP(SP1)!=1) {
1452 else if( (NbPoints==1)&&(NbPointsTotal==2) ) {
1453 if ( (SP1.CheckSameSP(Pt1))||(SP2.CheckSameSP(Pt1)) )
1455 else NbPointsTotal=3;
1457 else if( (NbPoints==2)&&(NbPointsTotal==0) ) {
1462 else if( (NbPoints==2)&&(NbPointsTotal==1) ) {//there is also Pt1 != Pt2
1463 if(SP1.CheckSameSP(Pt1)) {
1467 else if (SP1.CheckSameSP(Pt2)) {
1471 else NbPointsTotal=3;///case SP1!=Pt1 && SP1!=Pt2!
1473 else if( (NbPoints==2)&&(NbPointsTotal==2) ) {//there is also SP1!=SP2
1474 if( (SP1.CheckSameSP(Pt1))||(SP1.CheckSameSP(Pt2)) ) {
1475 if( (SP2.CheckSameSP(Pt1))||(SP2.CheckSameSP(Pt2)) )
1477 else NbPointsTotal=3;
1479 else NbPointsTotal=3;
1483 //=======================================================================
1484 //function : StartingPointsResearch
1486 //=======================================================================
1487 Standard_Integer IntPolyh_MaillageAffinage::StartingPointsResearch
1488 (const Standard_Integer T1,
1489 const Standard_Integer T2,
1490 IntPolyh_StartPoint &SP1,
1491 IntPolyh_StartPoint &SP2) const
1493 const IntPolyh_Point &P1=TPoints1[TTriangles1[T1].FirstPoint()];
1494 const IntPolyh_Point &P2=TPoints1[TTriangles1[T1].SecondPoint()];
1495 const IntPolyh_Point &P3=TPoints1[TTriangles1[T1].ThirdPoint()];
1496 const IntPolyh_Point &Q1=TPoints2[TTriangles2[T2].FirstPoint()];
1497 const IntPolyh_Point &Q2=TPoints2[TTriangles2[T2].SecondPoint()];
1498 const IntPolyh_Point &Q3=TPoints2[TTriangles2[T2].ThirdPoint()];
1501 /* The first triangle is (p1,p2,p3). The other is (q1,q2,q3).
1502 The sides are (e1,e2,e3) and (f1,f2,f3).
1503 The normals are n1 and m1*/
1505 const IntPolyh_Point e1=P2-P1;
1506 const IntPolyh_Point e2=P3-P2;
1507 const IntPolyh_Point e3=P1-P3;
1509 const IntPolyh_Point f1=Q2-Q1;
1510 const IntPolyh_Point f2=Q3-Q2;
1511 const IntPolyh_Point f3=Q1-Q3;
1514 IntPolyh_Point nn1,mm1;
1515 nn1.Cross(e1, e2); //normal of the first triangle
1516 mm1.Cross(f1, f2); //normal of the second triangle
1518 Standard_Real nn1modulus, mm1modulus;
1519 nn1modulus=sqrt(nn1.SquareModulus());
1520 mm1modulus=sqrt(mm1.SquareModulus());
1522 //-------------------------------------------------------
1523 ///calculation of intersection points between two triangles
1524 //-------------------------------------------------------
1525 Standard_Integer NbPoints=0;
1526 Standard_Integer NbPointsTotal=0;
1527 IntPolyh_StartPoint Pt1,Pt2;
1531 if(Abs(nn1modulus)<MyConfusionPrecision){//10.0e-20) {
1535 const IntPolyh_Point n1=nn1.Divide(nn1modulus);
1537 if(NbPointsTotal<2) {
1538 NbPoints=TriangleEdgeContact(1,1,P1,P2,P3,e1,e2,e3,Q1,Q2,f1,n1,Pt1,Pt2);
1539 TestNbPoints(1,NbPoints,NbPointsTotal,Pt1,Pt2,SP1,SP2);
1542 if(NbPointsTotal<2) {
1543 NbPoints=TriangleEdgeContact(1,2,P1,P2,P3,e1,e2,e3,Q2,Q3,f2,n1,Pt1,Pt2);
1544 TestNbPoints(1,NbPoints,NbPointsTotal,Pt1,Pt2,SP1,SP2);
1547 if(NbPointsTotal<2) {
1548 NbPoints=TriangleEdgeContact(1,3,P1,P2,P3,e1,e2,e3,Q3,Q1,f3,n1,Pt1,Pt2);
1549 TestNbPoints(1,NbPoints,NbPointsTotal,Pt1,Pt2,SP1,SP2);
1554 if(Abs(mm1modulus)<MyConfusionPrecision) {//10.0e-20){
1558 const IntPolyh_Point m1=mm1.Divide(mm1modulus);
1560 if(NbPointsTotal<2) {
1561 NbPoints=TriangleEdgeContact(2,1,Q1,Q2,Q3,f1,f2,f3,P1,P2,e1,m1,Pt1,Pt2);
1562 TestNbPoints(2,NbPoints,NbPointsTotal,Pt1,Pt2,SP1,SP2);
1565 if(NbPointsTotal<2) {
1566 NbPoints=TriangleEdgeContact(2,2,Q1,Q2,Q3,f1,f2,f3,P2,P3,e2,m1,Pt1,Pt2);
1567 TestNbPoints(2,NbPoints,NbPointsTotal,Pt1,Pt2,SP1,SP2);
1570 if(NbPointsTotal<2) {
1571 NbPoints=TriangleEdgeContact(2,3,Q1,Q2,Q3,f1,f2,f3,P3,P1,e3,m1,Pt1,Pt2);
1572 TestNbPoints(2,NbPoints,NbPointsTotal,Pt1,Pt2,SP1,SP2);
1576 /* if( (NbPointsTotal >1)&&( Abs(SP1.U1()-SP2.U1())<MyConfusionPrecision)
1577 &&(Abs(SP1.V1()-SP2.V1())<MyConfusionPrecision) )*/
1578 if( (NbPoints)&&(SP1.CheckSameSP(SP2)) )
1581 SP1.SetCoupleValue(T1,T2);
1582 SP2.SetCoupleValue(T1,T2);
1583 return (NbPointsTotal);
1585 //=======================================================================
1586 //function : StartingPointsResearch2
1587 //purpose : From two triangles compute intersection points.
1588 // If I found more than two intersection points
1589 // it means that those triangle are coplanar
1590 //=======================================================================
1591 Standard_Integer IntPolyh_MaillageAffinage::StartingPointsResearch2
1592 (const Standard_Integer T1,
1593 const Standard_Integer T2,
1594 IntPolyh_StartPoint &SP1,
1595 IntPolyh_StartPoint &SP2) const
1597 const IntPolyh_Triangle &Tri1=TTriangles1[T1];
1598 const IntPolyh_Triangle &Tri2=TTriangles2[T2];
1600 const IntPolyh_Point &P1=TPoints1[Tri1.FirstPoint()];
1601 const IntPolyh_Point &P2=TPoints1[Tri1.SecondPoint()];
1602 const IntPolyh_Point &P3=TPoints1[Tri1.ThirdPoint()];
1603 const IntPolyh_Point &Q1=TPoints2[Tri2.FirstPoint()];
1604 const IntPolyh_Point &Q2=TPoints2[Tri2.SecondPoint()];
1605 const IntPolyh_Point &Q3=TPoints2[Tri2.ThirdPoint()];
1609 /* The first triangle is (p1,p2,p3). The other is (q1,q2,q3).
1610 The sides are (e1,e2,e3) and (f1,f2,f3).
1611 The normals are n1 and m1*/
1613 const IntPolyh_Point e1=P2-P1;
1614 const IntPolyh_Point e2=P3-P2;
1615 const IntPolyh_Point e3=P1-P3;
1617 const IntPolyh_Point f1=Q2-Q1;
1618 const IntPolyh_Point f2=Q3-Q2;
1619 const IntPolyh_Point f3=Q1-Q3;
1622 IntPolyh_Point nn1,mm1;
1623 nn1.Cross(e1, e2); //normal to the first triangle
1624 mm1.Cross(f1, f2); //normal to the second triangle
1626 Standard_Real nn1modulus, mm1modulus;
1627 nn1modulus=sqrt(nn1.SquareModulus());
1628 mm1modulus=sqrt(mm1.SquareModulus());
1630 //-------------------------------------------------
1631 ///calculation of intersections points between triangles
1632 //-------------------------------------------------
1633 Standard_Integer NbPoints=0;
1634 Standard_Integer NbPointsTotal=0;
1638 if(Abs(nn1modulus)<MyConfusionPrecision) {//10.0e-20){
1642 const IntPolyh_Point n1=nn1.Divide(nn1modulus);
1644 if(NbPointsTotal<3) {
1645 IntPolyh_StartPoint Pt1,Pt2;
1646 NbPoints=TriangleEdgeContact2(1,1,Tri1,Tri2,P1,P2,P3,e1,e2,e3,Q1,Q2,f1,n1,Pt1,Pt2);
1647 TestNbPoints(1,NbPoints,NbPointsTotal,Pt1,Pt2,SP1,SP2);
1650 if(NbPointsTotal<3) {
1651 IntPolyh_StartPoint Pt1,Pt2;
1652 NbPoints=TriangleEdgeContact2(1,2,Tri1,Tri2,P1,P2,P3,e1,e2,e3,Q2,Q3,f2,n1,Pt1,Pt2);
1653 TestNbPoints(1,NbPoints,NbPointsTotal,Pt1,Pt2,SP1,SP2);
1656 if(NbPointsTotal<3) {
1657 IntPolyh_StartPoint Pt1,Pt2;
1658 NbPoints=TriangleEdgeContact2(1,3,Tri1,Tri2,P1,P2,P3,e1,e2,e3,Q3,Q1,f3,n1,Pt1,Pt2);
1659 TestNbPoints(1,NbPoints,NbPointsTotal,Pt1,Pt2,SP1,SP2);
1664 if(Abs(mm1modulus)<MyConfusionPrecision) {//10.0e-20){
1668 const IntPolyh_Point m1=mm1.Divide(mm1modulus);
1670 if(NbPointsTotal<3) {
1671 IntPolyh_StartPoint Pt1,Pt2;
1672 NbPoints=TriangleEdgeContact2(2,1,Tri1,Tri2,Q1,Q2,Q3,f1,f2,f3,P1,P2,e1,m1,Pt1,Pt2);
1673 TestNbPoints(2,NbPoints,NbPointsTotal,Pt1,Pt2,SP1,SP2);
1676 if(NbPointsTotal<3) {
1677 IntPolyh_StartPoint Pt1,Pt2;
1678 NbPoints=TriangleEdgeContact2(2,2,Tri1,Tri2,Q1,Q2,Q3,f1,f2,f3,P2,P3,e2,m1,Pt1,Pt2);
1679 TestNbPoints(2,NbPoints,NbPointsTotal,Pt1,Pt2,SP1,SP2);
1682 if(NbPointsTotal<3) {
1683 IntPolyh_StartPoint Pt1,Pt2;
1684 NbPoints=TriangleEdgeContact2(2,3,Tri1,Tri2,Q1,Q2,Q3,f1,f2,f3,P3,P1,e3,m1,Pt1,Pt2);
1685 TestNbPoints(2,NbPoints,NbPointsTotal,Pt1,Pt2,SP1,SP2);
1688 // no need already checked in TestNbPoints
1689 // if( (NbPointsTotal==2)&&(SP1.CheckSameSP(SP2)) ) {
1691 //SP1.SetCoupleValue(T1,T2);
1694 if(NbPointsTotal==2) {
1695 SP1.SetCoupleValue(T1,T2);
1696 SP2.SetCoupleValue(T1,T2);
1698 else if(NbPointsTotal==1)
1699 SP1.SetCoupleValue(T1,T2);
1700 else if(NbPointsTotal==3)
1701 SP1.SetCoupleValue(T1,T2);
1703 return (NbPointsTotal);
1705 //=======================================================================
1706 //function : NextStartingPointsResearch
1708 //=======================================================================
1709 Standard_Integer IntPolyh_MaillageAffinage::NextStartingPointsResearch
1710 (const Standard_Integer T1,
1711 const Standard_Integer T2,
1712 const IntPolyh_StartPoint &SPInit,
1713 IntPolyh_StartPoint &SPNext) const
1715 Standard_Integer NbPointsTotal=0;
1716 if( (T1<0)||(T2<0) ) NbPointsTotal=0;
1718 const IntPolyh_Point &P1=TPoints1[TTriangles1[T1].FirstPoint()];
1719 const IntPolyh_Point &P2=TPoints1[TTriangles1[T1].SecondPoint()];
1720 const IntPolyh_Point &P3=TPoints1[TTriangles1[T1].ThirdPoint()];
1721 const IntPolyh_Point &Q1=TPoints2[TTriangles2[T2].FirstPoint()];
1722 const IntPolyh_Point &Q2=TPoints2[TTriangles2[T2].SecondPoint()];
1723 const IntPolyh_Point &Q3=TPoints2[TTriangles2[T2].ThirdPoint()];
1725 /* The first triangle is (p1,p2,p3). The other is (q1,q2,q3).
1726 The sides are (e1,e2,e3) and (f1,f2,f3).
1727 The normals are n1 and m1*/
1729 const IntPolyh_Point e1=P2-P1;
1730 const IntPolyh_Point e2=P3-P2;
1731 const IntPolyh_Point e3=P1-P3;
1733 const IntPolyh_Point f1=Q2-Q1;
1734 const IntPolyh_Point f2=Q3-Q2;
1735 const IntPolyh_Point f3=Q1-Q3;
1737 IntPolyh_Point nn1,mm1;
1738 nn1.Cross(e1, e2); //normal to the first triangle
1739 mm1.Cross(f1, f2); //normal to the second triangle
1741 Standard_Real nn1modulus, mm1modulus;
1742 nn1modulus=sqrt(nn1.SquareModulus());
1743 mm1modulus=sqrt(mm1.SquareModulus());
1745 //-------------------------------------------------
1746 ///calculation of intersections points between triangles
1747 //-------------------------------------------------
1748 Standard_Integer NbPoints=0;
1749 IntPolyh_StartPoint SP1,SP2;
1752 if(Abs(nn1modulus)<MyConfusionPrecision) {//10.0e-20){
1756 const IntPolyh_Point n1=nn1.Divide(nn1modulus);
1758 if(NbPointsTotal<3) {
1759 IntPolyh_StartPoint Pt1,Pt2;
1760 NbPoints=TriangleEdgeContact(1,1,P1,P2,P3,e1,e2,e3,Q1,Q2,f1,n1,Pt1,Pt2);
1761 TestNbPoints(1,NbPoints,NbPointsTotal,Pt1,Pt2,SP1,SP2);
1764 if(NbPointsTotal<3) {
1765 IntPolyh_StartPoint Pt1,Pt2;
1766 NbPoints=TriangleEdgeContact(1,2,P1,P2,P3,e1,e2,e3,Q2,Q3,f2,n1,Pt1,Pt2);
1767 TestNbPoints(1,NbPoints,NbPointsTotal,Pt1,Pt2,SP1,SP2);
1770 if(NbPointsTotal<3) {
1771 IntPolyh_StartPoint Pt1,Pt2;
1772 NbPoints=TriangleEdgeContact(1,3,P1,P2,P3,e1,e2,e3,Q3,Q1,f3,n1,Pt1,Pt2);
1773 TestNbPoints(1,NbPoints,NbPointsTotal,Pt1,Pt2,SP1,SP2);
1778 if(Abs(mm1modulus)<MyConfusionPrecision) {//10.0e-20){
1782 const IntPolyh_Point m1=mm1.Divide(mm1modulus);
1784 if(NbPointsTotal<3) {
1785 IntPolyh_StartPoint Pt1,Pt2;
1786 NbPoints=TriangleEdgeContact(2,1,Q1,Q2,Q3,f1,f2,f3,P1,P2,e1,m1,Pt1,Pt2);
1787 TestNbPoints(2,NbPoints,NbPointsTotal,Pt1,Pt2,SP1,SP2);
1790 if(NbPointsTotal<3) {
1791 IntPolyh_StartPoint Pt1,Pt2;
1792 NbPoints=TriangleEdgeContact(2,2,Q1,Q2,Q3,f1,f2,f3,P2,P3,e2,m1,Pt1,Pt2);
1793 TestNbPoints(2,NbPoints,NbPointsTotal,Pt1,Pt2,SP1,SP2);
1796 if(NbPointsTotal<3) {
1797 IntPolyh_StartPoint Pt1,Pt2;
1798 NbPoints=TriangleEdgeContact(2,3,Q1,Q2,Q3,f1,f2,f3,P3,P1,e3,m1,Pt1,Pt2);
1799 TestNbPoints(2,NbPoints,NbPointsTotal,Pt1,Pt2,SP1,SP2);
1803 if (NbPointsTotal==1) {
1804 /* if( (Abs(SP1.U1()-SPInit.U1())<MyConfusionPrecision)
1805 &&(Abs(SP1.V1()-SPInit.V1())<MyConfusionPrecision) ) */
1806 if(SP1.CheckSameSP(SP2))
1814 // if ( (NbPointsTotal==2)&&( Abs(SP1.U1()-SPInit.U1())<MyConfusionPrecision)
1815 //&&( Abs(SP1.V1()-SPInit.V1())<MyConfusionPrecision) ) {
1816 if( (NbPointsTotal==2)&&(SP1.CheckSameSP(SPInit)) ) {
1817 NbPointsTotal=1;//SP1 et SPInit sont identiques
1820 // if( (NbPointsTotal==2)&&( Abs(SP2.U1()-SPInit.U1())<MyConfusionPrecision)
1821 //&&( Abs(SP2.V1()-SPInit.V1())<MyConfusionPrecision) ) {
1822 if( (NbPointsTotal==2)&&(SP2.CheckSameSP(SPInit)) ) {
1823 NbPointsTotal=1;//SP2 et SPInit sont identiques
1826 if(NbPointsTotal>1) {
1830 SPNext.SetCoupleValue(T1,T2);
1831 return (NbPointsTotal);
1833 //=======================================================================
1834 //function : NextStartingPointsResearch2
1835 //purpose : from two triangles and an intersection point I
1836 // seach the other point (if it exist).
1837 // This function is used by StartPointChain
1838 //=======================================================================
1839 Standard_Integer IntPolyh_MaillageAffinage::NextStartingPointsResearch2
1840 (const Standard_Integer T1,
1841 const Standard_Integer T2,
1842 const IntPolyh_StartPoint &SPInit,
1843 IntPolyh_StartPoint &SPNext) const
1845 Standard_Integer NbPointsTotal=0;
1846 Standard_Integer EdgeInit1=SPInit.E1();
1847 Standard_Integer EdgeInit2=SPInit.E2();
1848 if( (T1<0)||(T2<0) ) NbPointsTotal=0;
1851 const IntPolyh_Triangle &Tri1=TTriangles1[T1];
1852 const IntPolyh_Triangle &Tri2=TTriangles2[T2];
1854 const IntPolyh_Point &P1=TPoints1[Tri1.FirstPoint()];
1855 const IntPolyh_Point &P2=TPoints1[Tri1.SecondPoint()];
1856 const IntPolyh_Point &P3=TPoints1[Tri1.ThirdPoint()];
1857 const IntPolyh_Point &Q1=TPoints2[Tri2.FirstPoint()];
1858 const IntPolyh_Point &Q2=TPoints2[Tri2.SecondPoint()];
1859 const IntPolyh_Point &Q3=TPoints2[Tri2.ThirdPoint()];
1861 /* The first triangle is (p1,p2,p3). The other is (q1,q2,q3).
1862 The edges are (e1,e2,e3) and (f1,f2,f3).
1863 The normals are n1 and m1*/
1865 const IntPolyh_Point e1=P2-P1;
1866 const IntPolyh_Point e2=P3-P2;
1867 const IntPolyh_Point e3=P1-P3;
1869 const IntPolyh_Point f1=Q2-Q1;
1870 const IntPolyh_Point f2=Q3-Q2;
1871 const IntPolyh_Point f3=Q1-Q3;
1873 IntPolyh_Point nn1,mm1;
1874 nn1.Cross(e1, e2); //normal to the first triangle
1875 mm1.Cross(f1, f2); //normal to the second triangle
1877 Standard_Real nn1modulus, mm1modulus;
1878 nn1modulus=sqrt(nn1.SquareModulus());
1879 mm1modulus=sqrt(mm1.SquareModulus());
1881 //-------------------------------------------------
1882 ///calculation of intersections points between triangles
1883 //-------------------------------------------------
1885 Standard_Integer NbPoints=0;
1886 IntPolyh_StartPoint SP1,SP2;
1889 if(Abs(nn1modulus)<MyConfusionPrecision) {//10.0e-20){
1893 const IntPolyh_Point n1=nn1.Divide(nn1modulus);
1895 if( (NbPointsTotal<3)&&(EdgeInit2!=Tri2.FirstEdge()) ) {
1896 IntPolyh_StartPoint Pt1,Pt2;
1897 NbPoints=TriangleEdgeContact2(1,1,Tri1,Tri2,P1,P2,P3,e1,e2,e3,Q1,Q2,f1,n1,Pt1,Pt2);
1898 TestNbPoints(1,NbPoints,NbPointsTotal,Pt1,Pt2,SP1,SP2);
1901 if( (NbPointsTotal<3)&&(EdgeInit2!=Tri2.SecondEdge()) ) {
1902 IntPolyh_StartPoint Pt1,Pt2;
1903 NbPoints=TriangleEdgeContact2(1,2,Tri1,Tri2,P1,P2,P3,e1,e2,e3,Q2,Q3,f2,n1,Pt1,Pt2);
1904 TestNbPoints(1,NbPoints,NbPointsTotal,Pt1,Pt2,SP1,SP2);
1907 if( (NbPointsTotal<3)&&(EdgeInit2!=Tri2.ThirdEdge()) ) {
1908 IntPolyh_StartPoint Pt1,Pt2;
1909 NbPoints=TriangleEdgeContact2(1,3,Tri1,Tri2,P1,P2,P3,e1,e2,e3,Q3,Q1,f3,n1,Pt1,Pt2);
1910 TestNbPoints(1,NbPoints,NbPointsTotal,Pt1,Pt2,SP1,SP2);
1914 if(Abs(mm1modulus)<MyConfusionPrecision) {//10.0e-20){
1918 const IntPolyh_Point m1=mm1.Divide(mm1modulus);
1920 if( (NbPointsTotal<3)&&(EdgeInit1!=Tri1.FirstEdge()) ) {
1921 IntPolyh_StartPoint Pt1,Pt2;
1922 NbPoints=TriangleEdgeContact2(2,1,Tri1,Tri2,Q1,Q2,Q3,f1,f2,f3,P1,P2,e1,m1,Pt1,Pt2);
1923 TestNbPoints(2,NbPoints,NbPointsTotal,Pt1,Pt2,SP1,SP2);
1926 if( (NbPointsTotal<3)&&(EdgeInit1!=Tri1.SecondEdge()) ) {
1927 IntPolyh_StartPoint Pt1,Pt2;
1928 NbPoints=TriangleEdgeContact2(2,2,Tri1,Tri2,Q1,Q2,Q3,f1,f2,f3,P2,P3,e2,m1,Pt1,Pt2);
1929 TestNbPoints(2,NbPoints,NbPointsTotal,Pt1,Pt2,SP1,SP2);
1932 if( (NbPointsTotal<3)&&(EdgeInit1!=Tri1.ThirdEdge()) ) {
1933 IntPolyh_StartPoint Pt1,Pt2;
1934 NbPoints=TriangleEdgeContact2(2,3,Tri1,Tri2,Q1,Q2,Q3,f1,f2,f3,P3,P1,e3,m1,Pt1,Pt2);
1935 TestNbPoints(2,NbPoints,NbPointsTotal,Pt1,Pt2,SP1,SP2);
1939 if (NbPointsTotal==1) {
1940 if(SP1.CheckSameSP(SPInit))
1946 else if( (NbPointsTotal==2)&&(SP1.CheckSameSP(SPInit)) ) {
1947 NbPointsTotal=1;//SP1 et SPInit sont identiques
1950 else if( (NbPointsTotal==2)&&(SP2.CheckSameSP(SPInit)) ) {
1951 NbPointsTotal=1;//SP2 et SPInit sont identiques
1955 else if(NbPointsTotal>1) {
1959 SPNext.SetCoupleValue(T1,T2);
1960 return (NbPointsTotal);
1962 //=======================================================================
1963 //function : CalculPtsInterTriEdgeCoplanaires
1965 //=======================================================================
1966 void CalculPtsInterTriEdgeCoplanaires(const Standard_Integer TriSurfID,
1967 const IntPolyh_Point &NormaleTri,
1968 const IntPolyh_Point &PE1,
1969 const IntPolyh_Point &PE2,
1970 const IntPolyh_Point &Edge,
1971 const IntPolyh_Point &PT1,
1972 const IntPolyh_Point &PT2,
1973 const IntPolyh_Point &Cote,
1974 const Standard_Integer CoteIndex,
1975 IntPolyh_StartPoint &SP1,
1976 IntPolyh_StartPoint &SP2,
1977 Standard_Integer &NbPoints)
1979 IntPolyh_Point TestParalleles;
1980 TestParalleles.Cross(Edge,Cote);
1981 if(sqrt(TestParalleles.SquareModulus())<MyConfusionPrecision) {
1983 Per.Cross(NormaleTri,Cote);
1984 Standard_Real p1p = Per.Dot(PE1);
1985 Standard_Real p2p = Per.Dot(PE2);
1986 Standard_Real p0p = Per.Dot(PT1);
1987 if ( ( (p1p>=p0p)&&(p2p<=p0p) )||( (p1p<=p0p)&&(p2p>=p0p) ) ) {
1988 Standard_Real lambda=(p1p-p0p)/(p1p-p2p);
1989 if (lambda<-MyConfusionPrecision) {
1992 IntPolyh_Point PIE=PE1+Edge*lambda;
1993 Standard_Real alpha=RealLast();
1994 if(Cote.X()!=0) alpha=(PIE.X()-PT1.X())/Cote.X();
1995 else if (Cote.Y()!=0) alpha=(PIE.Y()-PT1.Y())/Cote.Y();
1996 else if (Cote.Z()!=0) alpha=(PIE.Z()-PT1.Z())/Cote.Z();
2000 if (alpha<-MyConfusionPrecision) {
2005 SP1.SetXYZ(PIE.X(),PIE.Y(),PIE.Z());
2007 SP1.SetUV1(PT1.U()+Cote.U()*alpha,PT1.V()+Cote.V()*alpha);
2008 SP1.SetUV2(PIE.U(),PIE.V());
2009 SP1.SetEdge1(CoteIndex);
2012 else if (TriSurfID==2) {
2013 SP1.SetUV1(PIE.U(),PIE.V());
2014 SP1.SetUV2(PT1.U()+Cote.U()*alpha,PT1.V()+Cote.V()*alpha);
2015 SP1.SetEdge2(CoteIndex);
2023 else if (NbPoints==1) {
2024 SP2.SetXYZ(PIE.X(),PIE.Y(),PIE.Z());
2026 SP2.SetUV1(PT1.U()+Cote.U()*alpha,PT1.V()+Cote.V()*alpha);
2027 SP2.SetUV2(PIE.U(),PIE.V());
2028 SP2.SetEdge1(CoteIndex);
2031 else if (TriSurfID==2) {
2032 SP2.SetUV1(PIE.U(),PIE.V());
2033 SP2.SetUV2(PT1.U()+Cote.U()*alpha,PT1.V()+Cote.V()*alpha);
2034 SP2.SetEdge2(CoteIndex);
2042 else if( (NbPoints>2)||(NbPoints<0) ) {
2048 else { //Cote et Edge paralleles, avec les rejections precedentes ils sont sur la meme droite
2049 //On projette les points sur cette droite
2050 Standard_Real pe1p= Cote.Dot(PE1);
2051 Standard_Real pe2p= Cote.Dot(PE2);
2052 Standard_Real pt1p= Cote.Dot(PT1);
2053 Standard_Real pt2p= Cote.Dot(PT2);
2055 IntPolyh_Point PEP1,PTP1,PEP2,PTP2;
2057 //PEP1 et PEP2 sont les points de contact entre le triangle et l'edge dans le repere UV de l'edge
2058 //PTP1 et PTP2 sont les correspondants respectifs a PEP1 et PEP2 dans le repere UV du triangle
2062 if ( (pt1p<pe1p)&&(pe1p<=pt2p) ) {
2064 PTP1=PT1+Cote*((pe1p-pt1p)/(pt2p-pt1p));
2068 PTP2=PT1+Cote*((pe2p-pt1p)/(pt2p-pt1p));
2072 PEP2=PE1+Edge*((pt1p-pe1p)/(pe2p-pe1p));
2077 else if( (pt2p<pe1p)&&(pe1p<=pt1p) ) {
2079 PTP1=PT1+Cote*((pt1p-pe1p)/(pt1p-pt2p));
2083 PTP2=PT1+Cote*((pe2p-pt1p)/(pt2p-pt1p));
2087 PEP2=PE1+Edge*((pt2p-pe1p)/(pe2p-pe1p));
2095 if ( (pt1p<pe2p)&&(pe2p<=pt2p) ) {
2097 PTP1=PT1+Cote*((pe2p-pt1p)/(pt2p-pt1p));
2101 PTP2=PT1+Cote*((pe1p-pt1p)/(pt2p-pt1p));
2105 PEP2=PE2+Edge*((pt1p-pe1p)/(pe2p-pe1p));
2110 else if( (pt2p<pe2p)&&(pe2p<=pt1p) ) {
2112 PTP1=PT1+Cote*((pt1p-pe2p)/(pt1p-pt2p));
2116 PTP2=PT1+Cote*((pe1p-pt1p)/(pt2p-pt1p));
2120 PEP2=PE1+Edge*((pt2p-pe1p)/(pe2p-pe1p));
2128 if (Abs(PEP1.U()-PEP2.U())<MyConfusionPrecision
2129 &&(Abs(PEP1.V()-PEP2.V())<MyConfusionPrecision) ) NbPoints=1;
2131 SP1.SetXYZ(PEP1.X(),PEP1.Y(),PEP1.Z());
2133 SP1.SetUV1(PTP1.U(),PTP1.V());
2134 SP1.SetUV2(PEP1.U(),PEP1.V());
2135 SP1.SetEdge1(CoteIndex);
2138 SP1.SetUV1(PEP1.U(),PTP1.V());
2139 SP1.SetUV2(PTP1.U(),PEP1.V());
2140 SP1.SetEdge2(CoteIndex);
2144 SP2.SetXYZ(PEP2.X(),PEP2.Y(),PEP2.Z());
2146 SP2.SetUV1(PTP2.U(),PTP2.V());
2147 SP2.SetUV2(PEP2.U(),PEP2.V());
2148 SP2.SetEdge1(CoteIndex);
2151 SP2.SetUV1(PEP2.U(),PTP2.V());
2152 SP2.SetUV2(PTP2.U(),PEP2.V());
2153 SP2.SetEdge2(CoteIndex);
2159 //=======================================================================
2160 //function : CalculPtsInterTriEdgeCoplanaires2
2162 //=======================================================================
2163 void CalculPtsInterTriEdgeCoplanaires2(const Standard_Integer TriSurfID,
2164 const IntPolyh_Point &NormaleTri,
2165 const IntPolyh_Triangle &Tri1,
2166 const IntPolyh_Triangle &Tri2,
2167 const IntPolyh_Point &PE1,
2168 const IntPolyh_Point &PE2,
2169 const IntPolyh_Point &Edge,
2170 const Standard_Integer EdgeIndex,
2171 const IntPolyh_Point &PT1,
2172 const IntPolyh_Point &PT2,
2173 const IntPolyh_Point &Cote,
2174 const Standard_Integer CoteIndex,
2175 IntPolyh_StartPoint &SP1,
2176 IntPolyh_StartPoint &SP2,
2177 Standard_Integer &NbPoints)
2179 IntPolyh_Point TestParalleles;
2180 TestParalleles.Cross(Edge,Cote);
2182 if(sqrt(TestParalleles.SquareModulus())>MyConfusionPrecision) {
2183 ///Edge and side are not parallel
2185 Per.Cross(NormaleTri,Cote);
2186 Standard_Real p1p = Per.Dot(PE1);
2187 Standard_Real p2p = Per.Dot(PE2);
2188 Standard_Real p0p = Per.Dot(PT1);
2189 ///The edge are PT1 are projected on the perpendicular of the side in the plane of the triangle
2190 if ( ( (p1p>=p0p)&&(p0p>=p2p) )||( (p1p<=p0p)&&(p0p<=p2p) ) ) {
2191 Standard_Real lambda=(p1p-p0p)/(p1p-p2p);
2192 if (lambda<-MyConfusionPrecision) {
2196 if (Abs(lambda)<MyConfusionPrecision)//lambda=0
2198 else if (Abs(lambda)>1.0-MyConfusionPrecision)//lambda=1
2201 PIE=PE1+Edge*lambda;
2203 Standard_Real alpha=RealLast();
2204 if(Cote.X()!=0) alpha=(PIE.X()-PT1.X())/Cote.X();
2205 else if (Cote.Y()!=0) alpha=(PIE.Y()-PT1.Y())/Cote.Y();
2206 else if (Cote.Z()!=0) alpha=(PIE.Z()-PT1.Z())/Cote.Z();
2211 if (alpha<-MyConfusionPrecision) {
2216 SP1.SetXYZ(PIE.X(),PIE.Y(),PIE.Z());
2218 if(Abs(alpha)<MyConfusionPrecision) {//alpha=0
2219 SP1.SetUV1(PT1.U(),PT1.V());
2220 SP1.SetUV1(PIE.U(),PIE.V());
2223 if(Abs(alpha)>1.0-MyConfusionPrecision) {//alpha=1
2224 SP1.SetUV1(PT2.U(),PT2.V());
2225 SP1.SetUV1(PIE.U(),PIE.V());
2229 SP1.SetUV1(PT1.U()+Cote.U()*alpha,PT1.V()+Cote.V()*alpha);
2230 SP1.SetUV2(PIE.U(),PIE.V());
2231 SP1.SetEdge1(Tri1.GetEdgeNumber(CoteIndex));
2232 if (Tri1.GetEdgeOrientation(CoteIndex)>0) SP1.SetLambda1(alpha);
2233 else SP1.SetLambda1(1.0-alpha);
2237 else if (TriSurfID==2) {
2238 if(Abs(alpha)<MyConfusionPrecision) {//alpha=0
2239 SP1.SetUV1(PT1.U(),PT1.V());
2240 SP1.SetUV1(PIE.U(),PIE.V());
2243 if(Abs(alpha)>1.0-MyConfusionPrecision) {//alpha=1
2244 SP1.SetUV1(PT2.U(),PT2.V());
2245 SP1.SetUV1(PIE.U(),PIE.V());
2249 SP1.SetUV1(PIE.U(),PIE.V());
2250 SP1.SetUV2(PT1.U()+Cote.U()*alpha,PT1.V()+Cote.V()*alpha);
2251 SP1.SetEdge2(Tri2.GetEdgeNumber(CoteIndex));
2252 if (Tri2.GetEdgeOrientation(CoteIndex)>0) SP1.SetLambda2(alpha);
2253 else SP1.SetLambda2(1.0-alpha);
2262 else if (NbPoints==1) {
2263 SP2.SetXYZ(PIE.X(),PIE.Y(),PIE.Z());
2265 if(Abs(alpha)<MyConfusionPrecision) {//alpha=0
2266 SP2.SetUV1(PT1.U(),PT1.V());
2267 SP2.SetUV1(PIE.U(),PIE.V());
2270 if(Abs(alpha)>1.0-MyConfusionPrecision) {//alpha=1
2271 SP2.SetUV1(PT2.U(),PT2.V());
2272 SP2.SetUV1(PIE.U(),PIE.V());
2276 SP2.SetUV1(PT1.U()+Cote.U()*alpha,PT1.V()+Cote.V()*alpha);
2277 SP2.SetUV2(PIE.U(),PIE.V());
2278 SP2.SetEdge1(Tri1.GetEdgeNumber(CoteIndex));
2279 if (Tri1.GetEdgeOrientation(CoteIndex)>0) SP2.SetLambda1(alpha);
2280 else SP2.SetLambda1(1.0-alpha);
2284 else if (TriSurfID==2) {
2285 if(Abs(alpha)<MyConfusionPrecision) {//alpha=0
2286 SP2.SetUV1(PT1.U(),PT1.V());
2287 SP2.SetUV1(PIE.U(),PIE.V());
2290 if(Abs(alpha)>1.0-MyConfusionPrecision) {//alpha=1
2291 SP2.SetUV1(PT2.U(),PT2.V());
2292 SP2.SetUV1(PIE.U(),PIE.V());
2296 SP2.SetUV1(PIE.U(),PIE.V());
2297 SP2.SetUV2(PT1.U()+Cote.U()*alpha,PT1.V()+Cote.V()*alpha);
2298 SP2.SetEdge2(Tri2.GetEdgeNumber(CoteIndex));
2299 if (Tri1.GetEdgeOrientation(CoteIndex)>0) SP2.SetLambda2(alpha);
2300 else SP2.SetLambda2(1.0-alpha);
2309 else if( (NbPoints>2)||(NbPoints<0) ) {
2316 //Side and Edge are parallel, with previous
2317 //rejections they are at the same side
2318 //The points are projected on that side
2319 Standard_Real pe1p= Cote.Dot(PE1);
2320 Standard_Real pe2p= Cote.Dot(PE2);
2321 Standard_Real pt1p= Cote.Dot(PT1);
2322 Standard_Real pt2p= Cote.Dot(PT2);
2323 Standard_Real lambda1=0., lambda2=0., alpha1=0., alpha2=0.;
2324 IntPolyh_Point PEP1,PTP1,PEP2,PTP2;
2327 if ( (pt1p<pe1p)&&(pe1p<=pt2p) ) {
2330 alpha1=((pe1p-pt1p)/(pt2p-pt1p));
2331 PTP1=PT1+Cote*alpha1;
2336 alpha2=((pe2p-pt1p)/(pt2p-pt1p));
2337 PTP2=PT1+Cote*alpha2;
2341 lambda2=((pt1p-pe1p)/(pe2p-pe1p));
2342 PEP2=PE1+Edge*lambda2;
2348 else if( (pt2p<pe1p)&&(pe1p<=pt1p) ) {
2351 alpha1=((pt1p-pe1p)/(pt1p-pt2p));
2352 PTP1=PT1+Cote*alpha1;
2357 alpha2=((pe2p-pt1p)/(pt2p-pt1p));
2358 PTP2=PT1+Cote*alpha2;
2362 lambda2=((pt2p-pe1p)/(pe2p-pe1p));
2363 PEP2=PE1+Edge*lambda2;
2372 if ( (pt1p<pe2p)&&(pe2p<=pt2p) ) {
2375 alpha1=((pe2p-pt1p)/(pt2p-pt1p));
2376 PTP1=PT1+Cote*alpha1;
2381 alpha2=((pe1p-pt1p)/(pt2p-pt1p));
2382 PTP2=PT1+Cote*alpha2;
2386 lambda2=((pt1p-pe1p)/(pe2p-pe1p));
2387 PEP2=PE2+Edge*lambda2;
2393 else if( (pt2p<pe2p)&&(pe2p<=pt1p) ) {
2396 alpha1=((pt1p-pe2p)/(pt1p-pt2p));
2397 PTP1=PT1+Cote*alpha1;
2402 alpha2=((pe1p-pt1p)/(pt2p-pt1p));
2403 PTP2=PT1+Cote*alpha2;
2407 lambda2=((pt2p-pe1p)/(pe2p-pe1p));
2408 PEP2=PE1+Edge*lambda2;
2417 SP1.SetXYZ(PEP1.X(),PEP1.Y(),PEP1.Z());
2418 if (TriSurfID==1) {///cote appartient a Tri1
2419 SP1.SetUV1(PTP1.U(),PTP1.V());
2420 SP1.SetUV2(PEP1.U(),PEP1.V());
2421 SP1.SetEdge1(Tri1.GetEdgeNumber(CoteIndex));
2423 if(Tri1.GetEdgeOrientation(CoteIndex)>0) SP1.SetLambda1(alpha1);
2424 else SP1.SetLambda1(1.0-alpha1);
2426 if(Tri2.GetEdgeOrientation(EdgeIndex)>0) SP1.SetLambda2(lambda1);
2427 else SP1.SetLambda2(1.0-lambda1);
2429 if (TriSurfID==2) {///cote appartient a Tri2
2430 SP1.SetUV1(PEP1.U(),PTP1.V());
2431 SP1.SetUV2(PTP1.U(),PEP1.V());
2432 SP1.SetEdge2(Tri2.GetEdgeNumber(CoteIndex));
2434 if(Tri2.GetEdgeOrientation(CoteIndex)>0) SP1.SetLambda1(alpha1);
2435 else SP1.SetLambda1(1.0-alpha1);
2437 if(Tri1.GetEdgeOrientation(EdgeIndex)>0) SP1.SetLambda2(lambda1);
2438 else SP1.SetLambda2(1.0-lambda1);
2441 //It is checked if PEP1!=PEP2
2442 if ( (NbPoints==2)&&(Abs(PEP1.U()-PEP2.U())<MyConfusionPrecision)
2443 &&(Abs(PEP1.V()-PEP2.V())<MyConfusionPrecision) ) NbPoints=1;
2445 SP2.SetXYZ(PEP2.X(),PEP2.Y(),PEP2.Z());
2447 SP2.SetUV1(PTP2.U(),PTP2.V());
2448 SP2.SetUV2(PEP2.U(),PEP2.V());
2449 SP2.SetEdge1(Tri1.GetEdgeNumber(CoteIndex));
2451 if(Tri1.GetEdgeOrientation(CoteIndex)>0) SP2.SetLambda1(alpha1);
2452 else SP2.SetLambda1(1.0-alpha1);
2454 if(Tri2.GetEdgeOrientation(EdgeIndex)>0) SP2.SetLambda2(lambda1);
2455 else SP2.SetLambda2(1.0-lambda1);
2458 SP2.SetUV1(PEP2.U(),PTP2.V());
2459 SP2.SetUV2(PTP2.U(),PEP2.V());
2460 SP2.SetEdge2(Tri2.GetEdgeNumber(CoteIndex));
2462 if(Tri1.GetEdgeOrientation(CoteIndex)>0) SP2.SetLambda1(alpha1);
2463 else SP2.SetLambda1(1.0-alpha1);
2465 if(Tri2.GetEdgeOrientation(EdgeIndex)>0) SP2.SetLambda2(lambda1);
2466 else SP2.SetLambda2(1.0-lambda1);
2471 //Filter if the point is placed on top, the edge is set to -1
2473 if(Abs(SP1.Lambda1())<MyConfusionPrecision)
2475 if(Abs(SP1.Lambda1()-1)<MyConfusionPrecision)
2477 if(Abs(SP1.Lambda2())<MyConfusionPrecision)
2479 if(Abs(SP1.Lambda2()-1)<MyConfusionPrecision)
2483 if(Abs(SP2.Lambda1())<MyConfusionPrecision)
2485 if(Abs(SP2.Lambda1()-1)<MyConfusionPrecision)
2487 if(Abs(SP2.Lambda2())<MyConfusionPrecision)
2489 if(Abs(SP2.Lambda2()-1)<MyConfusionPrecision)
2493 //=======================================================================
2494 //function : TriangleEdgeContact
2496 //=======================================================================
2497 Standard_Integer IntPolyh_MaillageAffinage::TriangleEdgeContact
2498 (const Standard_Integer TriSurfID,
2499 const Standard_Integer EdgeIndex,
2500 const IntPolyh_Point &PT1,
2501 const IntPolyh_Point &PT2,
2502 const IntPolyh_Point &PT3,
2503 const IntPolyh_Point &Cote12,
2504 const IntPolyh_Point &Cote23,
2505 const IntPolyh_Point &Cote31,
2506 const IntPolyh_Point &PE1,const IntPolyh_Point &PE2,
2507 const IntPolyh_Point &Edge,
2508 const IntPolyh_Point &NormaleT,
2509 IntPolyh_StartPoint &SP1,IntPolyh_StartPoint &SP2) const
2512 Standard_Real lambda =0.;
2513 Standard_Real alpha =0.;
2514 Standard_Real beta =0.;
2517 //The edge, on which the point is located, is known.
2519 SP1.SetEdge2(EdgeIndex);
2520 SP2.SetEdge2(EdgeIndex);
2522 else if (TriSurfID==2) {
2523 SP1.SetEdge1(EdgeIndex);
2524 SP2.SetEdge1(EdgeIndex);
2530 /**The edge is projected on the normal of the triangle if yes
2531 --> free intersection (point I)--> start point is found*/
2532 Standard_Integer NbPoints=0;
2534 if(NormaleT.SquareModulus()==0) {
2537 else if( (Cote12.SquareModulus()==0)
2538 ||(Cote23.SquareModulus()==0)
2539 ||(Cote31.SquareModulus()==0) ) {
2542 else if(Edge.SquareModulus()==0) {
2546 Standard_Real pe1 = NormaleT.Dot(PE1);
2547 Standard_Real pe2 = NormaleT.Dot(PE2);
2548 Standard_Real pt1 = NormaleT.Dot(PT1);
2550 // PE1I = lambda.Edge
2552 if( (Abs(pe1-pe2)<MyConfusionPrecision)&&(Abs(pe1-pt1)<MyConfusionPrecision) ) {
2553 //edge and triangle are coplanar (two contact points maximum)
2555 //The tops of the triangle are projected on the perpendicular of the edge
2557 IntPolyh_Point PerpEdge;
2558 PerpEdge.Cross(NormaleT,Edge);
2559 Standard_Real pp1 = PerpEdge.Dot(PT1);
2560 Standard_Real pp2 = PerpEdge.Dot(PT2);
2561 Standard_Real pp3 = PerpEdge.Dot(PT3);
2562 Standard_Real ppe1 = PerpEdge.Dot(PE1);
2564 if ( ( (pp1>ppe1)&&(pp2<=ppe1)&&(pp3<=ppe1) ) || ( (pp1<ppe1)&&(pp2>=ppe1)&&(pp3>=ppe1) ) ){
2565 //there are two sides (commun top PT1) that can cut the edge
2568 CalculPtsInterTriEdgeCoplanaires(TriSurfID,NormaleT,
2569 PE1,PE2,Edge,PT1,PT2,
2570 Cote12,1,SP1,SP2,NbPoints);
2572 if ( (NbPoints>1)&&(Abs(SP2.U1()-SP1.U1())<MyConfusionPrecision)
2573 &&(Abs(SP1.V1()-SP2.V1())<MyConfusionPrecision) ) NbPoints=1;
2577 CalculPtsInterTriEdgeCoplanaires(TriSurfID,NormaleT,
2578 PE1,PE2,Edge,PT3,PT1,
2579 Cote31,3,SP1,SP2,NbPoints);
2583 if ( (NbPoints>1)&&(Abs(SP1.U1()-SP2.U1())<MyConfusionPrecision)
2584 &&(Abs(SP1.V2()-SP2.V1())<MyConfusionPrecision) ) NbPoints=1;
2585 if (NbPoints>=2) return(NbPoints);
2587 else if ( ( ( (pp2>ppe1)&&(pp1<=ppe1)&&(pp3<=ppe1) ) || ( (pp2<ppe1)&&(pp1>=ppe1)&&(pp3>=ppe1) ) )
2589 //there are two sides (common top PT2) that can cut the edge
2592 CalculPtsInterTriEdgeCoplanaires(TriSurfID,NormaleT,
2593 PE1,PE2,Edge,PT1,PT2,
2594 Cote12,1,SP1,SP2,NbPoints);
2596 if ( (NbPoints>1)&&(Abs(SP2.U1()-SP1.U1())<MyConfusionPrecision)
2597 &&(Abs(SP1.V1()-SP2.V1())<MyConfusionPrecision) ) NbPoints=1;
2600 if(NbPoints<2) CalculPtsInterTriEdgeCoplanaires(TriSurfID,NormaleT,
2601 PE1,PE2,Edge,PT2,PT3,
2602 Cote23,2,SP1,SP2,NbPoints);
2604 if ( (NbPoints>1)&&(Abs(SP2.U1()-SP1.U1())<MyConfusionPrecision)
2605 &&(Abs(SP1.V1()-SP2.V1())<MyConfusionPrecision) ) NbPoints=1;
2606 if (NbPoints>=2) return(NbPoints);
2608 else if ( (( (pp3>ppe1)&&(pp1<=ppe1)&&(pp2<=ppe1) ) || ( (pp3<ppe1)&&(pp1>=ppe1)&&(pp2>=ppe1) ))
2610 //there are two sides (common top PT3) that can cut the edge
2613 CalculPtsInterTriEdgeCoplanaires(TriSurfID,NormaleT,
2614 PE1,PE2,Edge,PT3,PT1,Cote31,
2615 3,SP1,SP2,NbPoints);
2617 if ( (NbPoints>1)&&(Abs(SP2.U1()-SP1.U1())<MyConfusionPrecision)
2618 &&(Abs(SP1.V1()-SP2.V1())<MyConfusionPrecision) ) NbPoints=1;
2621 if (NbPoints<2) CalculPtsInterTriEdgeCoplanaires(TriSurfID,NormaleT,
2622 PE1,PE2,Edge,PT2,PT3,
2623 Cote23,2,SP1,SP2,NbPoints);
2625 if ( (NbPoints>1)&&(Abs(SP2.U1()-SP1.U1())<MyConfusionPrecision)
2626 &&(Abs(SP2.V1()-SP1.V1())<MyConfusionPrecision) ) NbPoints=1;
2627 if (NbPoints>=2) return(NbPoints);
2631 //------------------------------------------------------
2632 // edge and triangle NON COPLANAR (a contact point)
2633 //------------------------------------------------------
2634 else if( ( (pe1>=pt1)&&(pe2<=pt1) ) || ( (pe1<=pt1)&&(pe2>=pt1) ) ) {
2635 lambda=(pe1-pt1)/(pe1-pe2);
2637 if (lambda<-MyConfusionPrecision) {
2640 else if (Abs(lambda)<MyConfusionPrecision) {//lambda==0
2642 if(TriSurfID==1) SP1.SetEdge2(0);
2643 else SP1.SetEdge1(0);
2645 else if (Abs(lambda-1.0)<MyConfusionPrecision) {//lambda==1
2647 if(TriSurfID==1) SP1.SetEdge2(0);
2648 else SP1.SetEdge1(0);
2652 if(TriSurfID==1) SP1.SetEdge2(EdgeIndex);
2653 else SP1.SetEdge1(EdgeIndex);
2656 if(Abs(Cote23.X())>MyConfusionPrecision) {
2657 Standard_Real D=(Cote12.Y()-Cote12.X()*Cote23.Y()/Cote23.X());
2658 if(D!=0) alpha = (PI.Y()-PT1.Y()-(PI.X()-PT1.X())*Cote23.Y()/Cote23.X())/D;
2662 if ((alpha<-MyConfusionPrecision)||(alpha>(1.0+MyConfusionPrecision))) return(0);
2663 else beta = (PI.X()-PT1.X()-alpha*Cote12.X())/Cote23.X();
2666 else if (Abs(Cote12.X())>MyConfusionPrecision) { //On a Cote23.X()==0
2667 alpha = (PI.X()-PT1.X())/Cote12.X();
2669 if ((alpha<-MyConfusionPrecision)||(alpha>(1.0+MyConfusionPrecision))) return(0);
2671 else if (Abs(Cote23.Y())>MyConfusionPrecision) beta = (PI.Y()-PT1.Y()-alpha*Cote12.Y())/Cote23.Y();
2672 else if (Abs(Cote23.Z())>MyConfusionPrecision) beta = (PI.Z()-PT1.Z()-alpha*Cote12.Z())/Cote23.Z();
2678 else if (Abs(Cote23.Y())>MyConfusionPrecision) {
2679 //On a Cote23.X()==0 et Cote12.X()==0 ==> equation can't be used
2680 Standard_Real D=(Cote12.Z()-Cote12.Y()*Cote23.Z()/Cote23.Y());
2682 if(D!=0) alpha = (PI.Z()-PT1.Z()-(PI.Y()-PT1.Y())*Cote23.Z()/Cote23.Y())/D;
2687 if ( (alpha<-MyConfusionPrecision)||(alpha>(1.0+MyConfusionPrecision)) ) return(0);
2688 else beta = (PI.Y()-PT1.Y()-alpha*Cote12.Y())/Cote23.Y();
2691 else if (Abs(Cote12.Y())>MyConfusionPrecision) {
2692 //On a Cote23.X()==0, Cote12.X()==0 et Cote23.Y()==0
2693 alpha = (PI.Y()-PT1.Y())/Cote12.Y();
2695 if ((Abs(alpha)<MyConfusionPrecision)||(Abs(alpha-1.0)<MyConfusionPrecision)) return(0);
2697 else if (Abs(Cote23.Z())>MyConfusionPrecision) beta = (PI.Z()-PT1.Z()-alpha*Cote12.Z())/Cote23.Z();
2704 else { //two equations of three can't be used
2710 if( (beta<-MyConfusionPrecision)||(beta>(alpha+MyConfusionPrecision)) ) return(0);
2712 SP1.SetXYZ(PI.X(),PI.Y(),PI.Z());
2715 SP1.SetUV2(PI.U(),PI.V());
2716 SP1.SetUV1(PT1.U()+Cote12.U()*alpha+Cote23.U()*beta, PT1.V()+Cote12.V()*alpha+Cote23.V()*beta);
2718 if (beta<MyConfusionPrecision) {//beta==0 && alpha
2720 SP1.SetLambda1(alpha);
2722 if (Abs(beta-alpha)<MyConfusionPrecision) {//beta==alpha
2724 SP1.SetLambda1(1.0-alpha);
2726 if (Abs(alpha-1)<MyConfusionPrecision)//alpha==1
2728 if (alpha<MyConfusionPrecision) {//alpha=0 --> beta==0
2729 SP1.SetXYZ(PT1.X(),PT1.Y(),PT1.Z());
2730 SP1.SetUV1(PT1.U(),PT1.V());
2733 if ( (beta<MyConfusionPrecision)&&(Abs(1-alpha)<MyConfusionPrecision) ) {//beta==0 alpha==1
2734 SP1.SetXYZ(PT2.X(),PT2.Y(),PT2.Z());
2735 SP1.SetUV1(PT2.U(),PT2.V());
2738 if ( (Abs(beta-1)<MyConfusionPrecision)||(Abs(1-alpha)<MyConfusionPrecision) ) {//beta==1 alpha==1
2739 SP1.SetXYZ(PT3.X(),PT3.Y(),PT3.Z());
2740 SP1.SetUV1(PT3.U(),PT3.V());
2744 else if(TriSurfID==2) {
2745 SP1.SetUV1(PI.U(),PI.V());
2746 SP1.SetUV2(PT1.U()+Cote12.U()*alpha+Cote23.U()*beta, PT1.V()+Cote12.V()*alpha+Cote23.V()*beta);
2748 if (beta<MyConfusionPrecision) { //beta==0
2751 if (Abs(beta-alpha)<MyConfusionPrecision)//beta==alpha
2753 if (Abs(alpha-1)<MyConfusionPrecision)//alpha==1
2755 if (alpha<MyConfusionPrecision) {//alpha=0 --> beta==0
2756 SP1.SetXYZ(PT1.X(),PT1.Y(),PT1.Z());
2757 SP1.SetUV2(PT1.U(),PT1.V());
2760 if ( (beta<MyConfusionPrecision)&&(Abs(1-alpha)<MyConfusionPrecision) ) {//beta==0 alpha==1
2761 SP1.SetXYZ(PT2.X(),PT2.Y(),PT2.Z());
2762 SP1.SetUV2(PT2.U(),PT2.V());
2765 if ( (Abs(beta-1)<MyConfusionPrecision)||(Abs(1-alpha)<MyConfusionPrecision) ) {//beta==1 alpha==1
2766 SP1.SetXYZ(PT3.X(),PT3.Y(),PT3.Z());
2767 SP1.SetUV2(PT3.U(),PT3.V());
2781 //=======================================================================
2782 //function : TriangleEdgeContact2
2784 //=======================================================================
2785 Standard_Integer IntPolyh_MaillageAffinage::TriangleEdgeContact2
2786 (const Standard_Integer TriSurfID,
2787 const Standard_Integer EdgeIndex,
2788 const IntPolyh_Triangle &Tri1,
2789 const IntPolyh_Triangle &Tri2,
2790 const IntPolyh_Point &PT1,
2791 const IntPolyh_Point &PT2,
2792 const IntPolyh_Point &PT3,
2793 const IntPolyh_Point &Cote12,
2794 const IntPolyh_Point &Cote23,
2795 const IntPolyh_Point &Cote31,
2796 const IntPolyh_Point &PE1,const IntPolyh_Point &PE2,
2797 const IntPolyh_Point &Edge,
2798 const IntPolyh_Point &NormaleT,
2799 IntPolyh_StartPoint &SP1,IntPolyh_StartPoint &SP2) const
2802 Standard_Real lambda =0., alpha =0., beta =0.;
2804 //One of edges on which the point is located is known
2806 SP1.SetEdge2(Tri2.GetEdgeNumber(EdgeIndex));
2807 SP2.SetEdge2(Tri2.GetEdgeNumber(EdgeIndex));
2809 else if (TriSurfID==2) {
2810 SP1.SetEdge1(Tri1.GetEdgeNumber(EdgeIndex));
2811 SP2.SetEdge1(Tri1.GetEdgeNumber(EdgeIndex));
2817 /**The edge is projected on the normal in the triangle if yes
2818 --> a free intersection (point I)--> Start point is found */
2819 Standard_Integer NbPoints=0;
2820 if(NormaleT.SquareModulus()==0) {
2823 else if( (Cote12.SquareModulus()==0)
2824 ||(Cote23.SquareModulus()==0)
2825 ||(Cote31.SquareModulus()==0) ) {
2828 else if(Edge.SquareModulus()==0) {
2832 Standard_Real pe1 = NormaleT.Dot(PE1);
2833 Standard_Real pe2 = NormaleT.Dot(PE2);
2834 Standard_Real pt1 = NormaleT.Dot(PT1);
2836 // PE1I = lambda.Edge
2838 if( (Abs(pe1-pt1)<MyConfusionPrecision)&&(Abs(pe2-pt1)<MyConfusionPrecision)) {
2839 //edge and triangle are coplanar (two contact points at maximum)
2842 //the tops of the triangle are projected on the perpendicular to the edge
2843 IntPolyh_Point PerpEdge;
2844 PerpEdge.Cross(NormaleT,Edge);
2845 Standard_Real pp1 = PerpEdge.Dot(PT1);
2846 Standard_Real pp2 = PerpEdge.Dot(PT2);
2847 Standard_Real pp3 = PerpEdge.Dot(PT3);
2848 Standard_Real ppe1 = PerpEdge.Dot(PE1);
2851 if ( (Abs(pp1-pp2)<MyConfusionPrecision)&&(Abs(pp1-pp3)<MyConfusionPrecision) ) {
2855 if ( ( (pp1>=ppe1)&&(pp2<=ppe1)&&(pp3<=ppe1) ) || ( (pp1<=ppe1)&&(pp2>=ppe1)&&(pp3>=ppe1) ) ){
2856 //there are two sides (common top PT1) that can cut the edge
2859 CalculPtsInterTriEdgeCoplanaires2(TriSurfID,NormaleT,Tri1,Tri2,PE1,PE2,Edge,EdgeIndex,
2860 PT1,PT2,Cote12,1,SP1,SP2,NbPoints);
2862 if ( (NbPoints>1)&&(Abs(SP2.U1()-SP1.U1())<MyConfusionPrecision)
2863 &&(Abs(SP1.V1()-SP2.V1())<MyConfusionPrecision) ) NbPoints=1;
2866 if (NbPoints<2) CalculPtsInterTriEdgeCoplanaires2(TriSurfID,NormaleT,Tri1,Tri2,PE1,PE2,Edge,EdgeIndex,
2867 PT3,PT1,Cote31,3,SP1,SP2,NbPoints);
2870 if ( (NbPoints>1)&&(Abs(SP1.U1()-SP2.U1())<MyConfusionPrecision)
2871 &&(Abs(SP1.V2()-SP2.V1())<MyConfusionPrecision) ) NbPoints=1;
2872 if (NbPoints>=2) return(NbPoints);
2874 else if ( ( ( (pp2>=ppe1)&&(pp1<=ppe1)&&(pp3<=ppe1) ) || ( (pp2<=ppe1)&&(pp1>=ppe1)&&(pp3>=ppe1) ) )
2876 //there are two sides (common top PT2) that can cut the edge
2879 CalculPtsInterTriEdgeCoplanaires2(TriSurfID,NormaleT,Tri1,Tri2,PE1,PE2,Edge,EdgeIndex,
2880 PT1,PT2,Cote12,1,SP1,SP2,NbPoints);
2882 if ( (NbPoints>1)&&(Abs(SP2.U1()-SP1.U1())<MyConfusionPrecision)
2883 &&(Abs(SP1.V1()-SP2.V1())<MyConfusionPrecision) ) NbPoints=1;
2886 if(NbPoints<2) CalculPtsInterTriEdgeCoplanaires2(TriSurfID,NormaleT,Tri1,Tri2,PE1,PE2,Edge,EdgeIndex,
2887 PT2,PT3,Cote23,2,SP1,SP2,NbPoints);
2889 if ( (NbPoints>1)&&(Abs(SP2.U1()-SP1.U1())<MyConfusionPrecision)
2890 &&(Abs(SP1.V1()-SP2.V1())<MyConfusionPrecision) ) NbPoints=1;
2891 if (NbPoints>=2) return(NbPoints);
2893 else if ( (( (pp3>=ppe1)&&(pp1<=ppe1)&&(pp2<=ppe1) ) || ( (pp3<=ppe1)&&(pp1>=ppe1)&&(pp2>=ppe1) ))
2895 //there are two sides (common top PT3) that can cut the edge
2898 CalculPtsInterTriEdgeCoplanaires2(TriSurfID,NormaleT,Tri1,Tri2,PE1,PE2,Edge,EdgeIndex,
2899 PT3,PT1,Cote31,3,SP1,SP2,NbPoints);
2901 if ( (NbPoints>1)&&(Abs(SP2.U1()-SP1.U1())<MyConfusionPrecision)
2902 &&(Abs(SP1.V1()-SP2.V1())<MyConfusionPrecision) ) NbPoints=1;
2905 if (NbPoints<2) CalculPtsInterTriEdgeCoplanaires2(TriSurfID,NormaleT,Tri1,Tri2,PE1,PE2,Edge,EdgeIndex,
2906 PT2,PT3,Cote23,2,SP1,SP2,NbPoints);
2908 if ( (NbPoints>1)&&(Abs(SP2.U1()-SP1.U1())<MyConfusionPrecision)
2909 &&(Abs(SP2.V1()-SP1.V1())<MyConfusionPrecision) ) NbPoints=1;
2910 if (NbPoints>=2) return(NbPoints);
2914 //------------------------------------------------------
2915 // NON COPLANAR edge and triangle (a contact point)
2916 //------------------------------------------------------
2917 else if( ( (pe1>=pt1)&&(pt1>=pe2) ) || ( (pe1<=pt1)&&(pt1<=pe2) ) ) { //
2918 lambda=(pe1-pt1)/(pe1-pe2);
2920 if (lambda<-MyConfusionPrecision) {
2923 else if (Abs(lambda)<MyConfusionPrecision) {//lambda==0
2925 if(TriSurfID==1) SP1.SetEdge2(-1);
2926 else SP1.SetEdge1(-1);
2928 else if (Abs(lambda-1.0)<MyConfusionPrecision) {//lambda==1
2930 if(TriSurfID==1) SP1.SetEdge2(-1);
2931 else SP1.SetEdge1(-1);
2936 if(Tri2.GetEdgeOrientation(EdgeIndex)>0)
2937 SP1.SetLambda2(lambda);
2938 else SP1.SetLambda2(1.0-lambda);
2941 if(Tri1.GetEdgeOrientation(EdgeIndex)>0)
2942 SP1.SetLambda1(lambda);
2943 else SP1.SetLambda1(1.0-lambda);
2947 Standard_Real Cote23X=Cote23.X();
2948 Standard_Real D1=0.0;
2949 Standard_Real D3,D4;
2951 //Combination Eq1 Eq2
2952 if(Abs(Cote23X)>MyConfusionPrecision) {
2953 D1=Cote12.Y()-Cote12.X()*Cote23.Y()/Cote23X;
2955 if(Abs(D1)>MyConfusionPrecision) {
2956 alpha = ( PI.Y()-PT1.Y()-(PI.X()-PT1.X())*Cote23.Y()/Cote23X )/D1;
2958 ///It is checked if 1.0>=alpha>=0.0
2959 if ((alpha<-MyConfusionPrecision)||(alpha>(1.0+MyConfusionPrecision))) return(0);
2960 else beta = (PI.X()-PT1.X()-alpha*Cote12.X())/Cote23X;
2962 //Combination Eq1 and Eq2 with Cote23.X()==0
2963 else if ( (Abs(Cote12.X())>MyConfusionPrecision)
2964 &&(Abs(Cote23X)<MyConfusionPrecision) ) { //There is Cote23.X()==0
2965 alpha = (PI.X()-PT1.X())/Cote12.X();
2967 if ((alpha<-MyConfusionPrecision)||(alpha>(1.0+MyConfusionPrecision))) return(0);
2969 else if (Abs(Cote23.Y())>MyConfusionPrecision) beta = (PI.Y()-PT1.Y()-alpha*Cote12.Y())/Cote23.Y();
2970 else if (Abs(Cote23.Z())>MyConfusionPrecision) beta = (PI.Z()-PT1.Z()-alpha*Cote12.Z())/Cote23.Z();
2975 //Combination Eq1 and Eq3
2976 else if ( (Abs(Cote23.X())>MyConfusionPrecision)
2977 &&(Abs( D3= (Cote12.Z()-Cote12.X()*Cote23.Z()/Cote23.X()) ) > MyConfusionPrecision) ) {
2979 alpha = (PI.Z()-PT1.Z()-(PI.X()-PT1.X())*Cote23.Z()/Cote23.X())/D3;
2981 if ( (alpha<-MyConfusionPrecision)||(alpha>(1.0+MyConfusionPrecision)) ) return(0);
2982 else beta = (PI.X()-PT1.X()-alpha*Cote12.X())/Cote23.X();
2984 //Combination Eq2 and Eq3
2985 else if ( (Abs(Cote23.Y())>MyConfusionPrecision)
2986 &&(Abs( D4= (Cote12.Z()-Cote12.Y()*Cote23.Z()/Cote23.Y()) ) > MyConfusionPrecision) ) {
2988 alpha = (PI.Z()-PT1.Z()-(PI.Y()-PT1.Y())*Cote23.Z()/Cote23.Y())/D4;
2990 if ( (alpha<-MyConfusionPrecision)||(alpha>(1.0+MyConfusionPrecision)) ) return(0);
2991 else beta = (PI.Y()-PT1.Y()-alpha*Cote12.Y())/Cote23.Y();
2993 //Combination Eq2 and Eq3 with Cote23.Y()==0
2994 else if ( (Abs(Cote12.Y())>MyConfusionPrecision)
2995 && (Abs(Cote23.Y())<MyConfusionPrecision) ) {
2996 alpha = (PI.Y()-PT1.Y())/Cote12.Y();
2998 if ( (alpha<-MyConfusionPrecision)||(alpha>(1.0+MyConfusionPrecision)) ) return(0);
3000 else if (Abs(Cote23.Z())>MyConfusionPrecision) beta = (PI.Z()-PT1.Z()-alpha*Cote12.Z())/Cote23.Z();
3003 printf("\nCote PT2PT3 nul1\n");
3008 //Combination Eq1 and Eq3 with Cote23.Z()==0
3009 else if ( (Abs(Cote12.Z())>MyConfusionPrecision)
3010 && (Abs(Cote23.Z())<MyConfusionPrecision) ) {
3011 alpha = (PI.Z()-PT1.Z())/Cote12.Z();
3013 if ( (alpha<-MyConfusionPrecision)||(alpha>(1.0+MyConfusionPrecision)) ) return(0);
3015 else if (Abs(Cote23.X())>MyConfusionPrecision) beta = (PI.X()-PT1.X()-alpha*Cote12.X())/Cote23.X();
3022 else { //Particular case not processed ?
3028 if( (beta<-MyConfusionPrecision)||(beta>(alpha+MyConfusionPrecision)) ) return(0);
3030 SP1.SetXYZ(PI.X(),PI.Y(),PI.Z());
3033 SP1.SetUV2(PI.U(),PI.V());
3034 SP1.SetUV1(PT1.U()+Cote12.U()*alpha+Cote23.U()*beta, PT1.V()+Cote12.V()*alpha+Cote23.V()*beta);
3036 if (alpha<MyConfusionPrecision) {//alpha=0 --> beta==0
3037 SP1.SetXYZ(PT1.X(),PT1.Y(),PT1.Z());
3038 SP1.SetUV1(PT1.U(),PT1.V());
3041 else if ( (beta<MyConfusionPrecision)&&(Abs(1-alpha)<MyConfusionPrecision) ) {//beta==0 alpha==1
3042 SP1.SetXYZ(PT2.X(),PT2.Y(),PT2.Z());
3043 SP1.SetUV1(PT2.U(),PT2.V());
3046 else if ( (Abs(beta-1)<MyConfusionPrecision)&&(Abs(1-alpha)<MyConfusionPrecision) ) {//beta==1 alpha==1
3047 SP1.SetXYZ(PT3.X(),PT3.Y(),PT3.Z());
3048 SP1.SetUV1(PT3.U(),PT3.V());
3051 else if (beta<MyConfusionPrecision) {//beta==0
3052 SP1.SetEdge1(Tri1.GetEdgeNumber(1));
3053 if(Tri1.GetEdgeOrientation(1)>0)
3054 SP1.SetLambda1(alpha);
3055 else SP1.SetLambda1(1.0-alpha);
3057 else if (Abs(beta-alpha)<MyConfusionPrecision) {//beta==alpha
3058 SP1.SetEdge1(Tri1.GetEdgeNumber(3));
3059 if(Tri1.GetEdgeOrientation(3)>0)
3060 SP1.SetLambda1(1.0-alpha);
3061 else SP1.SetLambda1(alpha);
3063 else if (Abs(alpha-1)<MyConfusionPrecision) {//alpha==1
3064 SP1.SetEdge1(Tri1.GetEdgeNumber(2));
3065 if(Tri1.GetEdgeOrientation(2)>0)
3066 SP1.SetLambda1(beta);
3067 else SP1.SetLambda1(1.0-beta);
3070 else if(TriSurfID==2) {
3071 SP1.SetUV1(PI.U(),PI.V());
3072 SP1.SetUV2(PT1.U()+Cote12.U()*alpha+Cote23.U()*beta, PT1.V()+Cote12.V()*alpha+Cote23.V()*beta);
3074 if (alpha<MyConfusionPrecision) {//alpha=0 --> beta==0
3075 SP1.SetXYZ(PT1.X(),PT1.Y(),PT1.Z());
3076 SP1.SetUV2(PT1.U(),PT1.V());
3079 else if ( (beta<MyConfusionPrecision)&&(Abs(1-alpha)<MyConfusionPrecision) ) {//beta==0 alpha==1
3080 SP1.SetXYZ(PT2.X(),PT2.Y(),PT2.Z());
3081 SP1.SetUV2(PT2.U(),PT2.V());
3084 else if ( (Abs(beta-1)<MyConfusionPrecision)&&(Abs(1-alpha)<MyConfusionPrecision) ) {//beta==1 alpha==1
3085 SP1.SetXYZ(PT3.X(),PT3.Y(),PT3.Z());
3086 SP1.SetUV2(PT3.U(),PT3.V());
3089 else if (beta<MyConfusionPrecision) { //beta==0
3090 SP1.SetEdge2(Tri2.GetEdgeNumber(1));
3091 if(Tri2.GetEdgeOrientation(1)>0)
3092 SP1.SetLambda2(alpha);
3093 else SP1.SetLambda2(1.0-alpha);
3095 else if (Abs(beta-alpha)<MyConfusionPrecision) {//beta==alpha
3096 SP1.SetEdge2(Tri2.GetEdgeNumber(3));
3097 if(Tri2.GetEdgeOrientation(3)>0)
3098 SP1.SetLambda2(1.0-alpha);
3099 else SP1.SetLambda2(alpha);
3101 else if (Abs(alpha-1)<MyConfusionPrecision) {//alpha==1
3102 SP1.SetEdge2(Tri2.GetEdgeNumber(2));
3103 if(Tri2.GetEdgeOrientation(2)>0)
3104 SP1.SetLambda2(alpha);
3105 else SP1.SetLambda2(1.0-alpha);
3117 //=======================================================================
3118 //function : TriangleComparePSP
3119 //purpose : The same as TriangleCompare, plus compute the
3120 // StartPoints without chaining them.
3121 //=======================================================================
3122 Standard_Integer IntPolyh_MaillageAffinage::TriangleComparePSP ()
3124 Standard_Integer CpteurTab=0;
3125 Standard_Integer CpteurTabSP=0;
3126 Standard_Real CoupleAngle=-2.0;
3127 const Standard_Integer FinTT1 = TTriangles1.NbItems();
3128 const Standard_Integer FinTT2 = TTriangles2.NbItems();
3130 for(Standard_Integer i_S1=0; i_S1<FinTT1; i_S1++) {
3131 IntPolyh_Triangle &Triangle1 = TTriangles1[i_S1];
3132 if ((Triangle1.IndiceIntersectionPossible() == 0) ||
3133 (Triangle1.GetFleche() < 0.))
3135 for(Standard_Integer i_S2=0; i_S2<FinTT2; i_S2++){
3136 IntPolyh_Triangle &Triangle2 = TTriangles2[i_S2];
3137 if ((Triangle2.IndiceIntersectionPossible() != 0) &&
3138 (Triangle2.GetFleche() >= 0.)) {
3139 IntPolyh_StartPoint SP1, SP2;
3140 //If a triangle is dead or not in BSB, comparison is not possible
3142 Standard_Integer iDeg1, iDeg2, iDeg3, iDeg;
3144 const IntPolyh_Point& P1=TPoints1[Triangle1.FirstPoint()];
3145 const IntPolyh_Point& P2=TPoints1[Triangle1.SecondPoint()];
3146 const IntPolyh_Point& P3=TPoints1[Triangle1.ThirdPoint()];
3147 iDeg1=(P1.Degenerated()) ? 1 : 0;
3148 iDeg2=(P2.Degenerated()) ? 1 : 0;
3149 iDeg3=(P3.Degenerated()) ? 1 : 0;
3150 iDeg=iDeg1+iDeg2+iDeg3;
3155 const IntPolyh_Point& Q1=TPoints2[Triangle2.FirstPoint()];
3156 const IntPolyh_Point& Q2=TPoints2[Triangle2.SecondPoint()];
3157 const IntPolyh_Point& Q3=TPoints2[Triangle2.ThirdPoint()];
3158 iDeg1=(Q1.Degenerated()) ? 1 : 0;
3159 iDeg2=(Q2.Degenerated()) ? 1 : 0;
3160 iDeg3=(Q3.Degenerated()) ? 1 : 0;
3161 iDeg=iDeg1+iDeg2+iDeg3;
3166 if (TriContact(P1, P2, P3, Q1, Q2, Q3, CoupleAngle)) {
3167 Triangle1.SetIndiceIntersection(1);//The triangle is cut by another
3168 Triangle2.SetIndiceIntersection(1);
3170 Standard_Integer NbPoints;
3171 NbPoints=StartingPointsResearch(i_S1,i_S2,SP1, SP2);
3177 if ( (NbPoints>0)&&(NbPoints<3) ) {
3178 SP1.SetCoupleValue(i_S1,i_S2);
3179 TStartPoints[CpteurTabSP]=SP1;
3186 SP2.SetCoupleValue(i_S1,i_S2);
3187 TStartPoints[CpteurTabSP]=SP2;
3201 return(CpteurTabSP);
3203 //=======================================================================
3204 //function : TriangleCompare
3205 //purpose : Analyze each couple of triangles from the two --
3206 // array of triangles, to see if they are in
3207 // contact, and compute the incidence. Then put
3208 // couples in contact in the array of couples
3209 //=======================================================================
3210 Standard_Integer IntPolyh_MaillageAffinage::TriangleCompare ()
3212 Standard_Integer CpteurTab=0;
3214 const Standard_Integer FinTT1 = TTriangles1.NbItems();
3215 const Standard_Integer FinTT2 = TTriangles2.NbItems();
3217 Standard_Integer TTClimit = 200;
3218 Standard_Integer NbTTC = FinTT1 * FinTT2 / 10;
3219 if (NbTTC < TTClimit)
3221 TTrianglesContacts.Init(NbTTC);
3223 //TTrianglesContacts.Init(FinTT1 * FinTT2 / 10);
3225 Standard_Real CoupleAngle=-2.0;
3226 for(Standard_Integer i_S1=0; i_S1<FinTT1; i_S1++) {
3227 IntPolyh_Triangle &Triangle1 = TTriangles1[i_S1];
3228 if ((Triangle1.IndiceIntersectionPossible() == 0) ||
3229 (Triangle1.GetFleche() < 0.))
3231 for(Standard_Integer i_S2=0; i_S2<FinTT2; i_S2++){
3232 IntPolyh_Triangle &Triangle2 = TTriangles2[i_S2];
3233 if ((Triangle2.IndiceIntersectionPossible() != 0) &&
3234 (Triangle2.GetFleche() >= 0.)) {
3235 //If a triangle is dead or not in BSB, comparison is not possible
3236 Standard_Integer iDeg1, iDeg2, iDeg3, iDeg;
3238 const IntPolyh_Point& P1=TPoints1[Triangle1.FirstPoint()];
3239 const IntPolyh_Point& P2=TPoints1[Triangle1.SecondPoint()];
3240 const IntPolyh_Point& P3=TPoints1[Triangle1.ThirdPoint()];
3241 iDeg1=(P1.Degenerated()) ? 1 : 0;
3242 iDeg2=(P2.Degenerated()) ? 1 : 0;
3243 iDeg3=(P3.Degenerated()) ? 1 : 0;
3244 iDeg=iDeg1+iDeg2+iDeg3;
3249 const IntPolyh_Point& Q1=TPoints2[Triangle2.FirstPoint()];
3250 const IntPolyh_Point& Q2=TPoints2[Triangle2.SecondPoint()];
3251 const IntPolyh_Point& Q3=TPoints2[Triangle2.ThirdPoint()];
3252 iDeg1=(Q1.Degenerated()) ? 1 : 0;
3253 iDeg2=(Q2.Degenerated()) ? 1 : 0;
3254 iDeg3=(Q3.Degenerated()) ? 1 : 0;
3255 iDeg=iDeg1+iDeg2+iDeg3;
3260 if (TriContact(P1, P2, P3, Q1, Q2, Q3, CoupleAngle)) {
3261 if (CpteurTab >= NbTTC)
3263 TTrianglesContacts.SetNbItems(CpteurTab);
3266 TTrianglesContacts[CpteurTab].SetCoupleValue(i_S1, i_S2);
3267 TTrianglesContacts[CpteurTab].SetAngleValue(CoupleAngle);
3268 //test TTrianglesContacts[CpteurTab].Dump(CpteurTab);
3270 Triangle1.SetIndiceIntersection(1);//The triangle is cut by another
3271 Triangle2.SetIndiceIntersection(1);
3277 TTrianglesContacts.SetNbItems(CpteurTab);
3282 //=======================================================================
3283 //function : StartPointsCalcul
3284 //purpose : From the array of couples compute all the start
3285 // points and display them on the screen
3286 //=======================================================================
3287 void IntPolyh_MaillageAffinage::StartPointsCalcul() const
3289 const Standard_Integer FinTTC = TTrianglesContacts.NbItems();
3290 // printf("StartPointsCalcul() from IntPolyh_MaillageAffinage.cxx : StartPoints:\n");
3291 for(Standard_Integer ii=0; ii<FinTTC; ii++) {
3292 IntPolyh_StartPoint SP1,SP2;
3293 Standard_Integer T1,T2;
3294 T1=TTrianglesContacts[ii].FirstValue();
3295 T2=TTrianglesContacts[ii].SecondValue();
3296 StartingPointsResearch(T1,T2,SP1,SP2);
3297 if ( (SP1.E1()!=-1)&&(SP1.E2()!=-1) ) SP1.Dump(ii);
3298 if ( (SP2.E1()!=-1)&&(SP2.E2()!=-1) ) SP2.Dump(ii);
3301 //=======================================================================
3302 //function : CheckCoupleAndGetAngle
3304 //=======================================================================
3305 Standard_Integer CheckCoupleAndGetAngle(const Standard_Integer T1,
3306 const Standard_Integer T2,
3307 Standard_Real& Angle,
3308 IntPolyh_ArrayOfCouples &TTrianglesContacts)
3310 Standard_Integer Test=0;
3311 const Standard_Integer FinTTC = TTrianglesContacts.NbItems();
3312 for (Standard_Integer oioi=0; oioi<FinTTC; oioi++) {
3313 IntPolyh_Couple TestCouple = TTrianglesContacts[oioi];
3314 if ( (TestCouple.FirstValue()==T1)&&(TestCouple.AnalyseFlagValue()!=1) ) {
3315 if (TestCouple.SecondValue()==T2) {
3317 TTrianglesContacts[oioi].SetAnalyseFlag(1);
3318 Angle=TTrianglesContacts[oioi].AngleValue();
3325 //=======================================================================
3326 //function : CheckCoupleAndGetAngle2
3328 //=======================================================================
3329 Standard_Integer CheckCoupleAndGetAngle2(const Standard_Integer T1,
3330 const Standard_Integer T2,
3331 const Standard_Integer T11,
3332 const Standard_Integer T22,
3333 Standard_Integer &CT11,
3334 Standard_Integer &CT22,
3335 Standard_Real & Angle,
3336 IntPolyh_ArrayOfCouples &TTrianglesContacts)
3338 ///couple T1 T2 is found in the list
3339 ///T11 and T22 are two other triangles implied in the contact edge edge
3340 /// CT11 couple( T1,T22) and CT22 couple (T2,T11)
3341 /// these couples will be marked if there is a start point
3342 Standard_Integer Test1=0;
3343 Standard_Integer Test2=0;
3344 Standard_Integer Test3=0;
3345 const Standard_Integer FinTTC = TTrianglesContacts.NbItems();
3346 for (Standard_Integer oioi=0; oioi<FinTTC; oioi++) {
3347 IntPolyh_Couple TestCouple = TTrianglesContacts[oioi];
3348 if( (Test1==0)||(Test2==0)||(Test3==0) ) {
3349 if ( (TestCouple.FirstValue()==T1)&&(TestCouple.AnalyseFlagValue()!=1) ) {
3350 if (TestCouple.SecondValue()==T2) {
3352 TTrianglesContacts[oioi].SetAnalyseFlag(1);
3353 Angle=TTrianglesContacts[oioi].AngleValue();
3355 else if (TestCouple.SecondValue()==T22) {
3358 Angle=TTrianglesContacts[oioi].AngleValue();
3361 else if( (TestCouple.FirstValue()==T11)&&(TestCouple.AnalyseFlagValue()!=1) ) {
3362 if (TestCouple.SecondValue()==T2) {
3365 Angle=TTrianglesContacts[oioi].AngleValue();
3374 //=======================================================================
3375 //function : CheckNextStartPoint
3376 //purpose : it is checked if the point is not a top
3377 // then it is stored in one or several valid arrays with
3378 // the proper list number
3379 //=======================================================================
3380 Standard_Integer CheckNextStartPoint(IntPolyh_SectionLine & SectionLine,
3381 IntPolyh_ArrayOfTangentZones & TTangentZones,
3382 IntPolyh_StartPoint & SP,
3383 const Standard_Boolean Prepend)//=Standard_False)
3385 Standard_Integer Test=1;
3386 if( (SP.E1()==-1)||(SP.E2()==-1) ) {
3387 //The tops of triangle are analyzed
3388 //It is checked if they are not in the array TTangentZones
3389 Standard_Integer FinTTZ=TTangentZones.NbItems();
3390 for(Standard_Integer uiui=0; uiui<FinTTZ; uiui++) {
3391 IntPolyh_StartPoint TestSP=TTangentZones[uiui];
3392 if ( (Abs(SP.U1()-TestSP.U1())<MyConfusionPrecision)
3393 &&(Abs(SP.V1()-TestSP.V1())<MyConfusionPrecision) ) {
3394 if ( (Abs(SP.U2()-TestSP.U2())<MyConfusionPrecision)
3395 &&(Abs(SP.V2()-TestSP.V2())<MyConfusionPrecision) ) {
3396 Test=0;//SP is already in the list of tops
3401 if (Test) {//the top does not belong to the list of TangentZones
3402 SP.SetChainList(-1);
3403 TTangentZones[FinTTZ]=SP;
3404 TTangentZones.IncrementNbItems();
3405 Test=0;//the examined point is a top
3410 SectionLine.Prepend(SP);
3412 SectionLine[SectionLine.NbStartPoints()]=SP;
3413 SectionLine.IncrementNbStartPoints();
3417 //if the point is not a top Test=1
3418 //The chain is continued
3421 //=======================================================================
3422 //function : StartPointsChain
3423 //purpose : Loop on the array of couples. Compute StartPoints.
3424 // Try to chain the StartPoints into SectionLines or
3425 // put the point in the ArrayOfTangentZones if
3426 // chaining it, is not possible.
3427 //=======================================================================
3428 Standard_Integer IntPolyh_MaillageAffinage::StartPointsChain
3429 (IntPolyh_ArrayOfSectionLines& TSectionLines,
3430 IntPolyh_ArrayOfTangentZones& TTangentZones)
3432 //Loop on the array of couples filled in the function COMPARE()
3433 const Standard_Integer FinTTC = TTrianglesContacts.NbItems();
3435 //Array of tops of triangles
3436 for(Standard_Integer IndexA=0; IndexA<FinTTC; IndexA++) {
3437 //First couple of triangles.
3438 //It is checked if the couple of triangles has not been already examined.
3439 if(TTrianglesContacts[IndexA].AnalyseFlagValue()!=1) {
3441 Standard_Integer SectionLineIndex=TSectionLines.NbItems();
3442 // fill last section line if still empty (eap)
3443 if (SectionLineIndex > 0
3445 TSectionLines[SectionLineIndex-1].NbStartPoints() == 0)
3446 SectionLineIndex -= 1;
3448 TSectionLines.IncrementNbItems();
3450 IntPolyh_SectionLine & MySectionLine=TSectionLines[SectionLineIndex];
3451 if (MySectionLine.GetN() == 0) // eap
3452 MySectionLine.Init(10000);//Initialisation of array of StartPoint
3454 Standard_Integer NbPoints=-1;
3455 Standard_Integer T1I, T2I;
3456 T1I = TTrianglesContacts[IndexA].FirstValue();
3457 T2I = TTrianglesContacts[IndexA].SecondValue();
3459 // Start points for the current couple are found
3460 IntPolyh_StartPoint SP1, SP2;
3461 NbPoints=StartingPointsResearch2(T1I,T2I,SP1, SP2);//first calculation
3462 TTrianglesContacts[IndexA].SetAnalyseFlag(1);//the couple is marked
3464 if(NbPoints==1) {// particular case top/triangle or edge/edge
3465 //the start point is input in the array
3466 SP1.SetChainList(SectionLineIndex);
3467 SP1.SetAngle(TTrianglesContacts[IndexA].AngleValue());
3468 //it is checked if the point is not atop of the triangle
3469 if(CheckNextStartPoint(MySectionLine,TTangentZones,SP1)) {
3470 IntPolyh_StartPoint SPNext1;
3471 Standard_Integer TestSP1=0;
3474 IntPolyh_StartPoint SP11;//=SP1;
3475 if(SP1.E1()>=0) { //&&(SP1.E2()!=-1) already tested if the point is not a top
3476 Standard_Integer NextTriangle1=0;
3477 if (TEdges1[SP1.E1()].FirstTriangle()!=T1I) NextTriangle1=TEdges1[SP1.E1()].FirstTriangle();
3478 else NextTriangle1=TEdges1[SP1.E1()].SecondTriangle();
3480 Standard_Real Angle=-2.0;
3481 if (CheckCoupleAndGetAngle(NextTriangle1,T2I,Angle,TTrianglesContacts)) {
3482 //it is checked if the couple exists and is marked
3483 Standard_Integer NbPoints11=0;
3484 NbPoints11=NextStartingPointsResearch2(NextTriangle1,T2I,SP1,SP11);
3485 if (NbPoints11==1) {
3486 SP11.SetChainList(SectionLineIndex);
3487 SP11.SetAngle(Angle);
3489 if(CheckNextStartPoint(MySectionLine,TTangentZones,SP11)) {
3490 Standard_Integer EndChainList=1;
3491 while (EndChainList!=0) {
3492 TestSP1=GetNextChainStartPoint(SP11,SPNext1,MySectionLine,TTangentZones);
3494 SPNext1.SetChainList(SectionLineIndex);
3495 if(CheckNextStartPoint(MySectionLine,TTangentZones,SPNext1))
3497 else EndChainList=0;
3499 else EndChainList=0; //There is no next point
3505 if(NbPoints11>1) {//The point is input in the array TTangentZones
3506 TTangentZones[TTangentZones.NbItems()]=SP11;//default list number = -1
3507 TTangentZones.IncrementNbItems();
3515 else if (SP1.E2()<0){
3518 //chain of the other side
3519 IntPolyh_StartPoint SP12;//=SP1;
3520 if (SP1.E2()>=0) { //&&(SP1.E1()!=-1) already tested
3521 Standard_Integer NextTriangle2;
3522 if (TEdges2[SP1.E2()].FirstTriangle()!=T2I) NextTriangle2=TEdges2[SP1.E2()].FirstTriangle();
3523 else NextTriangle2=TEdges2[SP1.E2()].SecondTriangle();
3525 Standard_Real Angle=-2.0;
3526 if(CheckCoupleAndGetAngle(T1I,NextTriangle2,Angle,TTrianglesContacts)) {
3527 Standard_Integer NbPoints12=0;
3528 NbPoints12=NextStartingPointsResearch2(T1I,NextTriangle2,SP1, SP12);
3529 if (NbPoints12==1) {
3531 SP12.SetChainList(SectionLineIndex);
3532 SP12.SetAngle(Angle);
3533 Standard_Boolean Prepend = Standard_True; // eap
3535 if(CheckNextStartPoint(MySectionLine,TTangentZones,SP12, Prepend)) {
3536 Standard_Integer EndChainList=1;
3537 while (EndChainList!=0) {
3538 TestSP1=GetNextChainStartPoint(SP12,SPNext1,
3539 MySectionLine,TTangentZones,
3542 SPNext1.SetChainList(SectionLineIndex);
3543 if(CheckNextStartPoint(MySectionLine,TTangentZones,SPNext1,Prepend))
3545 else EndChainList=0;
3547 else EndChainList=0; //there is no next point
3552 if(NbPoints12>1) {//The points are input in the array TTangentZones
3553 TTangentZones[TTangentZones.NbItems()]=SP12;//default list number = -1
3554 TTangentZones.IncrementNbItems();
3563 else if(SP1.E1()<0){
3568 else if(NbPoints==2) {
3569 //the start points are input in the array
3570 IntPolyh_StartPoint SPNext2;
3571 Standard_Integer TestSP2=0;
3572 Standard_Integer EndChainList=1;
3574 SP1.SetChainList(SectionLineIndex);
3575 SP1.SetAngle(TTrianglesContacts[IndexA].AngleValue());
3576 if(CheckNextStartPoint(MySectionLine,TTangentZones,SP1)) {
3579 while (EndChainList!=0) {
3580 TestSP2=GetNextChainStartPoint(SP1,SPNext2,MySectionLine,TTangentZones);
3582 SPNext2.SetChainList(SectionLineIndex);
3583 if(CheckNextStartPoint(MySectionLine,TTangentZones,SPNext2))
3585 else EndChainList=0;
3587 else EndChainList=0; //there is no next point
3591 SP2.SetChainList(SectionLineIndex);
3592 SP2.SetAngle(TTrianglesContacts[IndexA].AngleValue());
3593 Standard_Boolean Prepend = Standard_True; // eap
3595 if(CheckNextStartPoint(MySectionLine,TTangentZones,SP2,Prepend)) {
3597 //chain of the other side
3599 while (EndChainList!=0) {
3600 TestSP2=GetNextChainStartPoint(SP2,SPNext2,
3601 MySectionLine,TTangentZones,
3604 SPNext2.SetChainList(SectionLineIndex);
3605 if(CheckNextStartPoint(MySectionLine,TTangentZones,SPNext2,Prepend))
3607 else EndChainList=0;
3609 else EndChainList=0; //there is no next point
3614 else if( (NbPoints>2)&&(NbPoints<7) ) {
3615 //More than two start points
3616 //the start point is input in the table
3617 SP1.SetChainList(SectionLineIndex);
3618 CheckNextStartPoint(MySectionLine,TTangentZones,SP1);
3629 //=======================================================================
3630 //function : GetNextChainStartPoint
3631 //purpose : Mainly used by StartPointsChain(), this function
3632 // try to compute the next StartPoint.
3633 // GetNextChainStartPoint is used only if it is known that there are 2 contact points
3634 //=======================================================================
3635 Standard_Integer IntPolyh_MaillageAffinage::GetNextChainStartPoint
3636 (const IntPolyh_StartPoint & SP,
3637 IntPolyh_StartPoint & SPNext,
3638 IntPolyh_SectionLine & MySectionLine,
3639 IntPolyh_ArrayOfTangentZones & TTangentZones,
3640 const Standard_Boolean Prepend)
3642 Standard_Integer NbPoints=0;
3643 if( (SP.E1()>=0)&&(SP.E2()==-2) ) {
3644 //case if the point is on edge of T1
3645 Standard_Integer NextTriangle1;
3646 if (TEdges1[SP.E1()].FirstTriangle()!=SP.T1()) NextTriangle1=TEdges1[SP.E1()].FirstTriangle();
3648 NextTriangle1=TEdges1[SP.E1()].SecondTriangle();
3649 //If is checked if two triangles intersect
3650 Standard_Real Angle= -2.0;
3651 if (CheckCoupleAndGetAngle(NextTriangle1,SP.T2(),Angle,TTrianglesContacts)) {
3652 NbPoints=NextStartingPointsResearch2(NextTriangle1,SP.T2(),SP,SPNext);
3655 CheckNextStartPoint(MySectionLine,TTangentZones,SPNext,Prepend);
3662 SPNext.SetAngle(Angle);
3664 else NbPoints=0;//this couple does not intersect
3666 else if( (SP.E1()==-2)&&(SP.E2()>=0) ) {
3667 //case if the point is on edge of T2
3668 Standard_Integer NextTriangle2;
3669 if (TEdges2[SP.E2()].FirstTriangle()!=SP.T2()) NextTriangle2=TEdges2[SP.E2()].FirstTriangle();
3671 NextTriangle2=TEdges2[SP.E2()].SecondTriangle();
3672 Standard_Real Angle= -2.0;
3673 if (CheckCoupleAndGetAngle(SP.T1(),NextTriangle2,Angle,TTrianglesContacts)) {
3674 NbPoints=NextStartingPointsResearch2(SP.T1(),NextTriangle2,SP,SPNext);
3677 CheckNextStartPoint(MySectionLine,TTangentZones,SPNext,Prepend);
3684 SPNext.SetAngle(Angle);
3688 else if( (SP.E1()==-2)&&(SP.E2()==-2) ) {
3689 ///no edge is touched or cut
3693 else if( (SP.E1()>=0)&&(SP.E2()>=0) ) {
3694 ///the point is located on two edges
3695 Standard_Integer NextTriangle1;
3696 Standard_Integer CpleT11=-1;
3697 Standard_Integer CpleT22=-1;
3698 if (TEdges1[SP.E1()].FirstTriangle()!=SP.T1()) NextTriangle1=TEdges1[SP.E1()].FirstTriangle();
3700 NextTriangle1=TEdges1[SP.E1()].SecondTriangle();
3701 Standard_Integer NextTriangle2;
3702 if (TEdges2[SP.E2()].FirstTriangle()!=SP.T2()) NextTriangle2=TEdges2[SP.E2()].FirstTriangle();
3704 NextTriangle2=TEdges2[SP.E2()].SecondTriangle();
3705 Standard_Real Angle= -2.0;
3706 if (CheckCoupleAndGetAngle2(NextTriangle1,NextTriangle2,
3707 SP.T1(),SP.T2(),CpleT11,CpleT22,
3708 Angle,TTrianglesContacts)) {
3709 NbPoints=NextStartingPointsResearch2(NextTriangle1,NextTriangle2,SP,SPNext);
3712 ///The new point is checked
3713 if(CheckNextStartPoint(MySectionLine,TTangentZones,SPNext,Prepend)>0) {
3722 SPNext.SetAngle(Angle);
3723 //The couples (Ti,Tj) (Ti',Tj') are marked
3724 if (CpleT11>=0) TTrianglesContacts[CpleT11].SetAnalyseFlag(1);
3728 if (CpleT22>=0) TTrianglesContacts[CpleT22].SetAnalyseFlag(1);
3736 else if( (SP.E1()==-1)||(SP.E2()==-1) ) {
3737 ///the points are tops of triangle
3738 ///the point is atored in an intermediary array
3742 //=======================================================================
3743 //function : GetArrayOfPoints
3745 //=======================================================================
3746 const IntPolyh_ArrayOfPoints& IntPolyh_MaillageAffinage::GetArrayOfPoints
3747 (const Standard_Integer SurfID)const
3753 //=======================================================================
3754 //function : GetArrayOfEdges
3756 //=======================================================================
3757 const IntPolyh_ArrayOfEdges& IntPolyh_MaillageAffinage::GetArrayOfEdges
3758 (const Standard_Integer SurfID)const
3764 //=======================================================================
3765 //function : GetArrayOfTriangles
3767 //=======================================================================
3768 const IntPolyh_ArrayOfTriangles&
3769 IntPolyh_MaillageAffinage::GetArrayOfTriangles
3770 (const Standard_Integer SurfID)const{
3772 return(TTriangles1);
3773 return(TTriangles2);
3776 //=======================================================================
3779 //=======================================================================
3780 Bnd_Box IntPolyh_MaillageAffinage::GetBox(const Standard_Integer SurfID) const
3787 //=======================================================================
3788 //function : GetBoxDraw
3790 //=======================================================================
3791 void IntPolyh_MaillageAffinage::GetBoxDraw(const Standard_Integer SurfID)const
3793 Standard_Real x0,y0,z0,x1,y1,z1;
3795 MyBox1.Get(x0,y0,z0,x1,y1,z1);
3798 MyBox2.Get(x0,y0,z0,x1,y1,z1);
3801 //=======================================================================
3802 //function : GetArrayOfCouples
3804 //=======================================================================
3805 IntPolyh_ArrayOfCouples &IntPolyh_MaillageAffinage::GetArrayOfCouples()
3807 return TTrianglesContacts;
3809 //=======================================================================
3810 //function : SetEnlargeZone
3812 //=======================================================================
3813 void IntPolyh_MaillageAffinage::SetEnlargeZone(Standard_Boolean& EnlargeZone)
3815 myEnlargeZone = EnlargeZone;
3817 //=======================================================================
3818 //function : GetEnlargeZone
3820 //=======================================================================
3821 Standard_Boolean IntPolyh_MaillageAffinage::GetEnlargeZone() const
3823 return myEnlargeZone;
3825 //=======================================================================
3826 //function : DegeneratedIndex
3828 //=======================================================================
3829 void DegeneratedIndex(const TColStd_Array1OfReal& aXpars,
3830 const Standard_Integer aNbX,
3831 const Handle(Adaptor3d_HSurface)& aS,
3832 const Standard_Integer aIsoDirection,
3833 Standard_Integer& aI1,
3834 Standard_Integer& aI2)
3837 Standard_Boolean bDegX1, bDegX2;
3838 Standard_Real aDegX1, aDegX2, aTol2, aX;
3842 aTol2=MyTolerance*MyTolerance;
3844 if (aIsoDirection==1){ // V=const
3845 bDegX1=IsDegenerated(aS, 1, aTol2, aDegX1);
3846 bDegX2=IsDegenerated(aS, 2, aTol2, aDegX2);
3848 else if (aIsoDirection==2){ // U=const
3849 bDegX1=IsDegenerated(aS, 3, aTol2, aDegX1);
3850 bDegX2=IsDegenerated(aS, 4, aTol2, aDegX2);
3856 if (!(bDegX1 || bDegX2)) {
3860 for(i=1; i<=aNbX; ++i) {
3863 if (fabs(aX-aDegX1) < MyTolerance) {
3868 if (fabs(aX-aDegX2) < MyTolerance) {
3874 //=======================================================================
3875 //function : IsDegenerated
3877 //=======================================================================
3878 Standard_Boolean IsDegenerated(const Handle(Adaptor3d_HSurface)& aS,
3879 const Standard_Integer aIndex,
3880 const Standard_Real aTol2,
3881 Standard_Real& aDegX)
3883 Standard_Boolean bRet;
3884 Standard_Integer i, aNbP;
3885 Standard_Real aU, dU, aU1, aU2, aV, dV, aV1, aV2, aD2;
3888 bRet=Standard_False;
3892 aU1=aS->FirstUParameter();
3893 aU2=aS->LastUParameter();
3894 aV1=aS->FirstVParameter();
3895 aV2=aS->LastVParameter();
3897 if (aIndex<3) { // V=const
3902 dU=(aU2-aU1)/(aNbP-1);
3904 aP1=aS->Value(aU, aV);
3905 for (i=1; i<aNbP; ++i) {
3910 aP2=aS->Value(aU, aV);
3911 aD2=aP1.SquareDistance(aP2);
3925 dV=(aV2-aV1)/(aNbP-1);
3927 aP1=aS->Value(aU, aV);
3928 for (i=1; i<aNbP; ++i) {
3933 aP2=aS->Value(aU, aV);
3934 aD2=aP1.SquareDistance(aP2);
3946 //=======================================================================
3947 //function : EnlargeZone
3949 //=======================================================================
3950 void EnlargeZone(const Handle(Adaptor3d_HSurface)& MaSurface,
3956 if(MaSurface->GetType() == GeomAbs_BSplineSurface ||
3957 MaSurface->GetType() == GeomAbs_BezierSurface) {
3958 if((!MaSurface->IsUClosed() && !MaSurface->IsUPeriodic()) &&
3959 (Abs(u0) < 1.e+100 && Abs(u1) < 1.e+100) ) {
3960 Standard_Real delta_u = 0.01*Abs(u1 - u0);
3964 if((!MaSurface->IsVClosed() && !MaSurface->IsVPeriodic()) &&
3965 (Abs(v0) < 1.e+100 && Abs(v1) < 1.e+100) ) {
3966 Standard_Real delta_v = 0.01*Abs(v1 - v0);