1 // Created on: 1999-03-05
2 // Created by: Fabrice SERVANT
3 // Copyright (c) 1999-1999 Matra Datavision
4 // Copyright (c) 1999-2012 OPEN CASCADE SAS
6 // The content of this file is subject to the Open CASCADE Technology Public
7 // License Version 6.5 (the "License"). You may not use the content of this file
8 // except in compliance with the License. Please obtain a copy of the License
9 // at http://www.opencascade.org and read it completely before using this file.
11 // The Initial Developer of the Original Code is Open CASCADE S.A.S., having its
12 // main offices at: 1, place des Freres Montgolfier, 78280 Guyancourt, France.
14 // The Original Code and all software distributed under the License is
15 // distributed on an "AS IS" basis, without warranty of any kind, and the
16 // Initial Developer hereby disclaims all such warranties, including without
17 // limitation, any warranties of merchantability, fitness for a particular
18 // purpose or non-infringement. Please see the License for the specific terms
19 // and conditions governing the rights and limitations under the License.
22 // modified by Edward AGAPOV (eap) Tue Jan 22 2002 (bug occ53)
23 // - improve SectionLine table management (avoid memory reallocation)
24 // - some protection against arrays overflow
26 // modified by Edward AGAPOV (eap) Thu Feb 14 2002 (occ139)
27 // - make Section Line parts rightly connected (prepend 2nd part to the 1st)
28 // - TriangleShape() for debugging purpose
30 // Modified by skv - Thu Sep 25 17:42:42 2003 OCC567
31 // modified by ofv Thu Apr 8 14:58:13 2004 fip
34 #include <IntPolyh_MaillageAffinage.ixx>
36 #include <Precision.hxx>
40 #include <TColStd_ListIteratorOfListOfInteger.hxx>
42 #include <Bnd_BoundSortBox.hxx>
43 #include <Bnd_HArray1OfBox.hxx>
45 #include <IntCurveSurface_ThePolyhedronOfHInter.hxx>
47 #include <IntPolyh_ArrayOfCouples.hxx>
48 #include <IntPolyh_Edge.hxx>
49 #include <IntPolyh_Couple.hxx>
51 static Standard_Real MyTolerance=10.0e-7;
52 static Standard_Real MyConfusionPrecision=10.0e-12;
53 static Standard_Real SquareMyConfusionPrecision=10.0e-24;
56 inline Standard_Real maxSR(const Standard_Real a,
57 const Standard_Real b,
58 const Standard_Real c);
61 inline Standard_Real minSR(const Standard_Real a,
62 const Standard_Real b,
63 const Standard_Real c);
65 Standard_Integer project6(const IntPolyh_Point &ax,
66 const IntPolyh_Point &p1,
67 const IntPolyh_Point &p2,
68 const IntPolyh_Point &p3,
69 const IntPolyh_Point &q1,
70 const IntPolyh_Point &q2,
71 const IntPolyh_Point &q3);
73 void TestNbPoints(const Standard_Integer ,
74 Standard_Integer &NbPoints,
75 Standard_Integer &NbPointsTotal,
76 const IntPolyh_StartPoint &Pt1,
77 const IntPolyh_StartPoint &Pt2,
78 IntPolyh_StartPoint &SP1,
79 IntPolyh_StartPoint &SP2);
81 void CalculPtsInterTriEdgeCoplanaires(const Standard_Integer TriSurfID,
82 const IntPolyh_Point &NormaleTri,
83 const IntPolyh_Point &PE1,
84 const IntPolyh_Point &PE2,
85 const IntPolyh_Point &Edge,
86 const IntPolyh_Point &PT1,
87 const IntPolyh_Point &PT2,
88 const IntPolyh_Point &Cote,
89 const Standard_Integer CoteIndex,
90 IntPolyh_StartPoint &SP1,
91 IntPolyh_StartPoint &SP2,
92 Standard_Integer &NbPoints);
94 void CalculPtsInterTriEdgeCoplanaires2(const Standard_Integer TriSurfID,
95 const IntPolyh_Point &NormaleTri,
96 const IntPolyh_Triangle &Tri1,
97 const IntPolyh_Triangle &Tri2,
98 const IntPolyh_Point &PE1,
99 const IntPolyh_Point &PE2,
100 const IntPolyh_Point &Edge,
101 const Standard_Integer EdgeIndex,
102 const IntPolyh_Point &PT1,
103 const IntPolyh_Point &PT2,
104 const IntPolyh_Point &Cote,
105 const Standard_Integer CoteIndex,
106 IntPolyh_StartPoint &SP1,
107 IntPolyh_StartPoint &SP2,
108 Standard_Integer &NbPoints);
110 Standard_Integer CheckCoupleAndGetAngle(const Standard_Integer T1,
111 const Standard_Integer T2,
112 Standard_Real& Angle,
113 IntPolyh_ArrayOfCouples &TTrianglesContacts);
115 Standard_Integer CheckCoupleAndGetAngle2(const Standard_Integer T1,
116 const Standard_Integer T2,
117 const Standard_Integer T11,
118 const Standard_Integer T22,
119 Standard_Integer &CT11,
120 Standard_Integer &CT22,
121 Standard_Real & Angle,
122 IntPolyh_ArrayOfCouples &TTrianglesContacts);
124 Standard_Integer CheckNextStartPoint(IntPolyh_SectionLine & SectionLine,
125 IntPolyh_ArrayOfTangentZones & TTangentZones,
126 IntPolyh_StartPoint & SP,
127 const Standard_Boolean Prepend=Standard_False);
130 Standard_Boolean IsDegenerated(const Handle(Adaptor3d_HSurface)& aS,
131 const Standard_Integer aIndex,
132 const Standard_Real aTol2,
133 Standard_Real& aDegX);
135 void DegeneratedIndex(const TColStd_Array1OfReal& Xpars,
136 const Standard_Integer aNbX,
137 const Handle(Adaptor3d_HSurface)& aS,
138 const Standard_Integer aIsoDirection,
139 Standard_Integer& aI1,
140 Standard_Integer& aI2);
142 void EnlargeZone(const Handle(Adaptor3d_HSurface)& MaSurface,
148 //=======================================================================
149 //function : IntPolyh_MaillageAffinage
151 //=======================================================================
152 IntPolyh_MaillageAffinage::IntPolyh_MaillageAffinage
153 (const Handle(Adaptor3d_HSurface)& Surface1,
154 const Handle(Adaptor3d_HSurface)& Surface2,
155 const Standard_Integer )
157 MaSurface1(Surface1),
158 MaSurface2(Surface2),
169 myEnlargeZone(Standard_False)
172 //=======================================================================
173 //function : IntPolyh_MaillageAffinage
175 //=======================================================================
176 IntPolyh_MaillageAffinage::IntPolyh_MaillageAffinage
177 (const Handle(Adaptor3d_HSurface)& Surface1,
178 const Standard_Integer NbSU1,
179 const Standard_Integer NbSV1,
180 const Handle(Adaptor3d_HSurface)& Surface2,
181 const Standard_Integer NbSU2,
182 const Standard_Integer NbSV2,
183 const Standard_Integer )
185 MaSurface1(Surface1),
186 MaSurface2(Surface2),
197 myEnlargeZone(Standard_False)
200 //=======================================================================
201 //function : FillArrayOfPnt
202 //purpose : Compute points on one surface and fill an array of points
203 //=======================================================================
204 void IntPolyh_MaillageAffinage::FillArrayOfPnt
205 (const Standard_Integer SurfID)
207 Standard_Integer NbSamplesU, NbSamplesV, i, aNbSamplesU1, aNbSamplesV1;
208 Standard_Real u0, u1, v0, v1, aU, aV, dU, dV;
210 const Handle(Adaptor3d_HSurface&) MaSurface=(SurfID==1)? MaSurface1 : MaSurface2;
211 NbSamplesU=(SurfID==1)? NbSamplesU1:NbSamplesU2;
212 NbSamplesV=(SurfID==1)? NbSamplesV1:NbSamplesV2;
214 u0 = (MaSurface)->FirstUParameter();
215 u1 = (MaSurface)->LastUParameter();
216 v0 = (MaSurface)->FirstVParameter();
217 v1 = (MaSurface)->LastVParameter();
220 EnlargeZone(MaSurface, u0, u1, v0, v1);
223 TColStd_Array1OfReal aUpars(1, NbSamplesU);
224 TColStd_Array1OfReal aVpars(1, NbSamplesV);
226 aNbSamplesU1=NbSamplesU-1;
227 aNbSamplesV1=NbSamplesV-1;
229 dU=(u1-u0)/Standard_Real(aNbSamplesU1);
230 dV=(v1-v0)/Standard_Real(aNbSamplesV1);
232 for (i=0; i<NbSamplesU; ++i) {
234 if (i==aNbSamplesU1) {
237 aUpars.SetValue(i+1, aU);
240 for (i=0; i<NbSamplesV; ++i) {
242 if (i==aNbSamplesV1) {
245 aVpars.SetValue(i+1, aV);
248 FillArrayOfPnt(SurfID, aUpars, aVpars);
250 //=======================================================================
251 //function : FillArrayOfPnt
252 //purpose : Compute points on one surface and fill an array of points
253 // FILL AN ARRAY OF POINTS
254 //=======================================================================
255 void IntPolyh_MaillageAffinage::FillArrayOfPnt
256 (const Standard_Integer SurfID,
257 const Standard_Boolean isShiftFwd)
259 Standard_Integer NbSamplesU, NbSamplesV, i, aNbSamplesU1, aNbSamplesV1;
260 Standard_Real u0, u1, v0, v1, aU, aV, dU, dV;
261 const Handle(Adaptor3d_HSurface)& MaSurface=(SurfID==1)? MaSurface1 : MaSurface2;
262 NbSamplesU=(SurfID==1)? NbSamplesU1:NbSamplesU2;
263 NbSamplesV=(SurfID==1)? NbSamplesV1:NbSamplesV2;
265 u0 = (MaSurface)->FirstUParameter();
266 u1 = (MaSurface)->LastUParameter();
267 v0 = (MaSurface)->FirstVParameter();
268 v1 = (MaSurface)->LastVParameter();
271 EnlargeZone(MaSurface, u0, u1, v0, v1);
274 TColStd_Array1OfReal aUpars(1, NbSamplesU);
275 TColStd_Array1OfReal aVpars(1, NbSamplesV);
277 aNbSamplesU1=NbSamplesU-1;
278 aNbSamplesV1=NbSamplesV-1;
280 dU=(u1-u0)/Standard_Real(aNbSamplesU1);
281 dV=(v1-v0)/Standard_Real(aNbSamplesV1);
283 for (i=0; i<NbSamplesU; ++i) {
285 if (i==aNbSamplesU1) {
288 aUpars.SetValue(i+1, aU);
291 for (i=0; i<NbSamplesV; ++i) {
293 if (i==aNbSamplesV1) {
296 aVpars.SetValue(i+1, aV);
299 FillArrayOfPnt(SurfID, isShiftFwd, aUpars, aVpars);
301 //=======================================================================
302 //function : FillArrayOfPnt
303 //purpose : Compute points on one surface and fill an array of points
304 //=======================================================================
305 void IntPolyh_MaillageAffinage::FillArrayOfPnt
306 (const Standard_Integer SurfID,
307 const TColStd_Array1OfReal& Upars,
308 const TColStd_Array1OfReal& Vpars)
310 Standard_Boolean bDegI, bDeg;
311 Standard_Integer aNbU, aNbV, iCnt, i, j;
312 Standard_Integer aID1, aID2, aJD1, aJD2;
313 Standard_Real aTol, aU, aV, aX, aY, aZ;
316 aNbU=(SurfID==1)? NbSamplesU1 : NbSamplesU2;
317 aNbV=(SurfID==1)? NbSamplesV1 : NbSamplesV2;
318 Bnd_Box& aBox = (SurfID==1) ? MyBox1 : MyBox2;
319 Handle(Adaptor3d_HSurface)& aS=(SurfID==1)? MaSurface1:MaSurface2;
320 IntPolyh_ArrayOfPoints &TPoints=(SurfID==1)? TPoints1:TPoints2;
326 DegeneratedIndex(Vpars, aNbV, aS, 1, aJD1, aJD2);
327 if (!(aJD1 || aJD2)) {
328 DegeneratedIndex(Upars, aNbU, aS, 2, aID1, aID2);
331 TPoints.Init(aNbU*aNbV);
333 for(i=1; i<=aNbU; ++i){
334 bDegI=(aID1==i || aID2==i);
336 for(j=1; j<=aNbV; ++j){
338 aP=aS->Value(aU, aV);
339 aP.Coord(aX, aY, aZ);
340 IntPolyh_Point& aIP=TPoints[iCnt];
341 aIP.Set(aX, aY, aZ, aU, aV);
343 bDeg=bDegI || (aJD1==j || aJD2==j);
345 aIP.SetDegenerated(bDeg);
352 TPoints.SetNbItems(iCnt);
354 IntCurveSurface_ThePolyhedronOfHInter polyhedron(aS, Upars, Vpars);
356 aTol=polyhedron.DeflectionOverEstimation();
359 Standard_Real a1,a2,a3,b1,b2,b3;
361 aBox.Get(a1,a2,a3,b1,b2,b3);
362 aBox.Update(a1-aTol,a2-aTol,a3-aTol,b1+aTol,b2+aTol,b3+aTol);
363 aBox.Enlarge(MyTolerance);
366 //=======================================================================
367 //function : FillArrayOfPnt
368 //purpose : Compute points on one surface and fill an array of points
369 // REMPLISSAGE DU TABLEAU DE POINTS
370 //=======================================================================
371 void IntPolyh_MaillageAffinage::FillArrayOfPnt
372 (const Standard_Integer SurfID,
373 const Standard_Boolean isShiftFwd,
374 const TColStd_Array1OfReal& Upars,
375 const TColStd_Array1OfReal& Vpars)
377 Standard_Boolean bDegI, bDeg;
378 Standard_Integer aNbU, aNbV, iCnt, i, j;
379 Standard_Integer aID1, aID2, aJD1, aJD2;
380 Standard_Real Tol, resol, u0, v0, u1, v1, aU, aV, aMag;
381 Standard_Real aX, aY, aZ;
383 gp_Vec aDU, aDV, aNorm;
385 aNbU=(SurfID==1)? NbSamplesU1:NbSamplesU2;
386 aNbV=(SurfID==1)? NbSamplesV1:NbSamplesV2;
387 Bnd_Box& aBox = (SurfID==1) ? MyBox1 : MyBox2;
388 Handle(Adaptor3d_HSurface) aS=(SurfID==1)? MaSurface1:MaSurface2;
389 IntPolyh_ArrayOfPoints &TPoints=(SurfID==1)? TPoints1:TPoints2;
391 resol = gp::Resolution();
397 IntCurveSurface_ThePolyhedronOfHInter polyhedron(aS, Upars, Vpars);
398 Tol=polyhedron.DeflectionOverEstimation();
403 DegeneratedIndex(Vpars, aNbV, aS, 1, aJD1, aJD2);
404 if (!(aJD1 || aJD2)) {
405 DegeneratedIndex(Upars, aNbU, aS, 2, aID1, aID2);
408 TPoints.Init(aNbU*aNbV);
410 for(i=1; i<=aNbU; ++i){
411 bDegI=(aID1==i || aID2==i);
413 for(j=1; j<=aNbV; ++j){
415 aS->D1(aU, aV, aP, aDU, aDV);
417 aNorm = aDU.Crossed(aDV);
418 aMag = aNorm.Magnitude();
421 aNorm.Multiply(Tol*1.5);
427 aP.Translate(aNorm.Reversed());
431 IntPolyh_Point& aIP=TPoints[iCnt];
432 aP.Coord(aX, aY, aZ);
433 aIP.Set(aX, aY, aZ, aU, aV);
435 bDeg=bDegI || (aJD1==j || aJD2==j);
437 aIP.SetDegenerated(bDeg);
444 TPoints.SetNbItems(iCnt);
448 Standard_Real a1,a2,a3,b1,b2,b3;
450 aBox.Get(a1,a2,a3,b1,b2,b3);
451 aBox.Update(a1-Tol,a2-Tol,a3-Tol,b1+Tol,b2+Tol,b3+Tol);
452 aBox.Enlarge(MyTolerance);
454 //=======================================================================
455 //function : CommonBox
456 //purpose : Compute the common box witch is the intersection
457 // of the two bounding boxes, and mark the points of
458 // the two surfaces that are inside.
459 // REJECTION BOUNDING BOXES
460 // DETERMINATION OF THE COMMON BOX
461 //=======================================================================
462 void IntPolyh_MaillageAffinage::CommonBox (const Bnd_Box &,
471 Standard_Real x10,y10,z10,x11,y11,z11;
472 Standard_Real x20,y20,z20,x21,y21,z21;
474 MyBox1.Get(x10,y10,z10,x11,y11,z11);
475 MyBox2.Get(x20,y20,z20,x21,y21,z21);
483 if((x10>x21)||(x20>x11)||(y10>y21)||(y20>y11)||(z10>z21)||(z20>z11)) {
486 if(x11<=x21) XMax=x11; else { if(x21<=x11) XMax=x21;}
487 if(x20<=x10) XMin=x10; else { if(x10<=x20) XMin=x20;}
488 if(y11<=y21) YMax=y11; else { if(y21<=y11) YMax=y21;}
489 if(y20<=y10) YMin=y10; else { if(y10<=y20) YMin=y20;}
490 if(z11<=z21) ZMax=z11; else { if(z21<=z11) ZMax=z21;}
491 if(z20<=z10) ZMin=z10; else { if(z10<=z20) ZMin=z20;}
493 if(((XMin==XMax)&&(!(YMin==YMax)&&!(ZMin==ZMax)))
494 ||((YMin==YMax)&&(!(XMin==XMax)&&!(ZMin==ZMax)))//ou exclusif ??
495 ||((ZMin==ZMax)&&(!(XMin==XMax)&&!(YMin==YMax)))) {
503 //extension of the box
504 if( (X==0)&&(Y!=0) ) X=Y*0.1;
505 else if( (X==0)&&(Z!=0) ) X=Z*0.1;
508 if( (Y==0)&&(X!=0) ) Y=X*0.1;
509 else if( (Y==0)&&(Z!=0) ) Y=Z*0.1;
512 if( (Z==0)&&(X!=0) ) Z=X*0.1;
513 else if( (Z==0)&&(Y!=0) ) Z=Y*0.1;
517 if( (X==0)&&(Y==0)&&(Z==0) ) {
525 //Marking of points included in the common
526 const Standard_Integer FinTP1 = TPoints1.NbItems();
527 // for(Standard_Integer i=0; i<FinTP1; i++) {
529 for( i=0; i<FinTP1; i++) {
530 IntPolyh_Point & Pt1 = TPoints1[i];
557 Pt1.SetPartOfCommon(r);
560 const Standard_Integer FinTP2 = TPoints2.NbItems();
561 for(Standard_Integer ii=0; ii<FinTP2; ii++) {
562 IntPolyh_Point & Pt2 = TPoints2[ii];
590 Pt2.SetPartOfCommon(rr);
593 //=======================================================================
594 //function : FillArrayOfEdges
595 //purpose : Compute edges from the array of points
596 // FILL THE ARRAY OF EDGES
597 //=======================================================================
598 void IntPolyh_MaillageAffinage::FillArrayOfEdges
599 (const Standard_Integer SurfID)
602 IntPolyh_ArrayOfEdges &TEdges=(SurfID==1)? TEdges1:TEdges2;
603 Standard_Integer NbSamplesU=(SurfID==1)? NbSamplesU1:NbSamplesU2;
604 Standard_Integer NbSamplesV=(SurfID==1)? NbSamplesV1:NbSamplesV2;
606 //NbEdges = 3 + 3*(NbSamplesV-2) + 3*(NbSamplesU-2) +
607 // + 3*(NbSamplesU-2)*(NbSamplesV-2) + (NbSamplesV-1) + (NbSamplesU-1);
608 //NbSamplesU and NbSamples cannot be less than 2, so
609 Standard_Integer NbEdges = 3*NbSamplesU*NbSamplesV - 2*(NbSamplesU+NbSamplesV) + 1;
610 TEdges.Init(NbEdges);
612 Standard_Integer CpteurTabEdges=0;
615 TEdges[CpteurTabEdges].SetFirstPoint(0); // U V
616 TEdges[CpteurTabEdges].SetSecondPoint(1); // U V+1
617 // TEdges[CpteurTabEdges].SetFirstTriangle(-1);
618 TEdges[CpteurTabEdges].SetSecondTriangle(0);
621 TEdges[CpteurTabEdges].SetFirstPoint(0); // U V
622 TEdges[CpteurTabEdges].SetSecondPoint(NbSamplesV); // U+1 V
623 TEdges[CpteurTabEdges].SetFirstTriangle(0);
624 TEdges[CpteurTabEdges].SetSecondTriangle(1);
627 TEdges[CpteurTabEdges].SetFirstPoint(0); // U V
628 TEdges[CpteurTabEdges].SetSecondPoint(NbSamplesV+1); // U+1 V+1
629 TEdges[CpteurTabEdges].SetFirstTriangle(1);
630 // TEdges[CpteurTabEdges].SetSecondTriangle(-1);
634 Standard_Integer PntInit=1;
635 Standard_Integer BoucleMeshV;
636 for(BoucleMeshV=1; BoucleMeshV<NbSamplesV-1;BoucleMeshV++){
637 TEdges[CpteurTabEdges].SetFirstPoint(PntInit); // U V
638 TEdges[CpteurTabEdges].SetSecondPoint(PntInit+1); // U V+1
639 // TEdges[CpteurTabEdges].SetFirstTriangle(-1);
640 TEdges[CpteurTabEdges].SetSecondTriangle(BoucleMeshV*2);
643 TEdges[CpteurTabEdges].SetFirstPoint(PntInit); // U V
644 TEdges[CpteurTabEdges].SetSecondPoint(PntInit+NbSamplesV+1); // U+1 V+1
645 TEdges[CpteurTabEdges].SetFirstTriangle(BoucleMeshV*2);
646 TEdges[CpteurTabEdges].SetSecondTriangle(BoucleMeshV*2+1);
649 TEdges[CpteurTabEdges].SetFirstPoint(PntInit); // U V
650 TEdges[CpteurTabEdges].SetSecondPoint(PntInit+NbSamplesV); // U+1 V
651 TEdges[CpteurTabEdges].SetFirstTriangle(BoucleMeshV*2+1);
652 TEdges[CpteurTabEdges].SetSecondTriangle(BoucleMeshV*2-2);
659 for(BoucleMeshV=1; BoucleMeshV<NbSamplesU-1;BoucleMeshV++){
660 TEdges[CpteurTabEdges].SetFirstPoint(PntInit); // U V
661 TEdges[CpteurTabEdges].SetSecondPoint(PntInit+1); // U V+1
662 TEdges[CpteurTabEdges].SetFirstTriangle((BoucleMeshV-1)*(NbSamplesV-1)*2+1);
663 TEdges[CpteurTabEdges].SetSecondTriangle(BoucleMeshV*(NbSamplesV-1)*2);
666 TEdges[CpteurTabEdges].SetFirstPoint(PntInit); // U V
667 TEdges[CpteurTabEdges].SetSecondPoint(PntInit+NbSamplesV+1); // U+1 V+1
668 TEdges[CpteurTabEdges].SetFirstTriangle(BoucleMeshV*(NbSamplesV-1)*2);
669 TEdges[CpteurTabEdges].SetSecondTriangle(BoucleMeshV*(NbSamplesV-1)*2+1);
672 TEdges[CpteurTabEdges].SetFirstPoint(PntInit); // U V
673 TEdges[CpteurTabEdges].SetSecondPoint(PntInit+NbSamplesV); // U+1 V
674 TEdges[CpteurTabEdges].SetFirstTriangle(BoucleMeshV*(NbSamplesV-1)*2+1);
675 // TEdges[CpteurTabEdges].SetSecondTriangle(-1);
680 PntInit=NbSamplesV+1;
681 //To provide recursion I associate a point with three edges
682 for(Standard_Integer BoucleMeshU=1; BoucleMeshU<NbSamplesU-1; BoucleMeshU++){
683 for(Standard_Integer BoucleMeshV=1; BoucleMeshV<NbSamplesV-1;BoucleMeshV++){
684 TEdges[CpteurTabEdges].SetFirstPoint(PntInit); // U V
685 TEdges[CpteurTabEdges].SetSecondPoint(PntInit+1); // U V+1
686 TEdges[CpteurTabEdges].SetFirstTriangle((NbSamplesV-1)*2*(BoucleMeshU-1)+BoucleMeshV*2+1);
687 TEdges[CpteurTabEdges].SetSecondTriangle((NbSamplesV-1)*2*BoucleMeshU+BoucleMeshV*2);
690 TEdges[CpteurTabEdges].SetFirstPoint(PntInit); // U V
691 TEdges[CpteurTabEdges].SetSecondPoint(PntInit+NbSamplesV+1); // U+1 V+1
692 TEdges[CpteurTabEdges].SetFirstTriangle((NbSamplesV-1)*2*BoucleMeshU+BoucleMeshV*2);
693 TEdges[CpteurTabEdges].SetSecondTriangle((NbSamplesV-1)*2*BoucleMeshU+BoucleMeshV*2+1);
696 TEdges[CpteurTabEdges].SetFirstPoint(PntInit); // U V
697 TEdges[CpteurTabEdges].SetSecondPoint(PntInit+NbSamplesV); // U+1 V
698 TEdges[CpteurTabEdges].SetFirstTriangle((NbSamplesV-1)*2*BoucleMeshU+BoucleMeshV*2+1);
699 TEdges[CpteurTabEdges].SetSecondTriangle((NbSamplesV-1)*2*BoucleMeshU+BoucleMeshV*2-2);
701 PntInit++;//Pass to the next point
703 PntInit++;//Pass the last point of the column
704 PntInit++;//Pass the first point of the next column
708 PntInit=(NbSamplesU-1)*NbSamplesV; //point U=u1 V=0
709 for(BoucleMeshV=0; BoucleMeshV<NbSamplesV-1; BoucleMeshV++){
710 TEdges[CpteurTabEdges].SetFirstPoint(PntInit); //U=u1 V
711 TEdges[CpteurTabEdges].SetSecondPoint(PntInit+1); //U=u1 V+1
712 TEdges[CpteurTabEdges].SetFirstTriangle((NbSamplesU-2)*(NbSamplesV-1)*2+BoucleMeshV*2+1);
713 // TEdges[CpteurTabEdges].SetSecondTriangle(-1);
719 for(BoucleMeshV=0; BoucleMeshV<NbSamplesU-1;BoucleMeshV++){
720 TEdges[CpteurTabEdges].SetFirstPoint(NbSamplesV-1+BoucleMeshV*NbSamplesV); // U V=v1
721 TEdges[CpteurTabEdges].SetSecondPoint(NbSamplesV-1+(BoucleMeshV+1)*NbSamplesV); //U+1 V=v1
722 // TEdges[CpteurTabEdges].SetFirstTriangle(-1);
723 TEdges[CpteurTabEdges].SetSecondTriangle(BoucleMeshV*2*(NbSamplesV-1)+(NbSamplesV-2)*2);
726 TEdges.SetNbItems(CpteurTabEdges);
729 //=======================================================================
730 //function : FillArrayOfTriangles
731 //purpose : Compute triangles from the array of points, and --
732 // mark the triangles that use marked points by the
733 // CommonBox function.
734 // FILL THE ARRAY OF TRIANGLES
735 //=======================================================================
736 void IntPolyh_MaillageAffinage::FillArrayOfTriangles
737 (const Standard_Integer SurfID)
739 Standard_Integer CpteurTabTriangles=0;
740 Standard_Integer PntInit=0;
742 IntPolyh_ArrayOfPoints &TPoints=(SurfID==1)? TPoints1:TPoints2;
743 IntPolyh_ArrayOfTriangles &TTriangles=(SurfID==1)? TTriangles1:TTriangles2;
744 Standard_Integer NbSamplesU=(SurfID==1)? NbSamplesU1:NbSamplesU2;
745 Standard_Integer NbSamplesV=(SurfID==1)? NbSamplesV1:NbSamplesV2;
747 TTriangles.Init(2*(NbSamplesU-1)*(NbSamplesV-1));
748 //To provide recursion, I associate a point with two triangles
749 for(Standard_Integer BoucleMeshU=0; BoucleMeshU<NbSamplesU-1; BoucleMeshU++){
750 for(Standard_Integer BoucleMeshV=0; BoucleMeshV<NbSamplesV-1;BoucleMeshV++){
753 TTriangles[CpteurTabTriangles].SetFirstPoint(PntInit); // U V
754 TTriangles[CpteurTabTriangles].SetSecondPoint(PntInit+1); // U V+1
755 TTriangles[CpteurTabTriangles].SetThirdPoint(PntInit+NbSamplesV+1); // U+1 V+1
757 // IF ITS EDGE CONTACTS WITH THE COMMON BOX IP REMAINS = A 1
758 if( ( (TPoints[PntInit].PartOfCommon()) & (TPoints[PntInit+1].PartOfCommon()) )
759 &&( (TPoints[PntInit+1].PartOfCommon()) & (TPoints[PntInit+NbSamplesV+1].PartOfCommon()))
760 &&( (TPoints[PntInit+NbSamplesV+1].PartOfCommon()) & (TPoints[PntInit].PartOfCommon())) )
762 TTriangles[CpteurTabTriangles].SetIndiceIntersectionPossible(0);
764 CpteurTabTriangles++;
767 TTriangles[CpteurTabTriangles].SetFirstPoint(PntInit); // U V
768 TTriangles[CpteurTabTriangles].SetSecondPoint(PntInit+NbSamplesV+1); // U+1 V+1
769 TTriangles[CpteurTabTriangles].SetThirdPoint(PntInit+NbSamplesV); // U+1 V
772 if( ( (TPoints[PntInit].PartOfCommon()) & (TPoints[PntInit+NbSamplesV+1].PartOfCommon()) )
773 &&( (TPoints[PntInit+NbSamplesV+1].PartOfCommon()) & (TPoints[PntInit+NbSamplesV].PartOfCommon()))
774 &&( (TPoints[PntInit+NbSamplesV].PartOfCommon()) & (TPoints[PntInit].PartOfCommon())) )
775 TTriangles[CpteurTabTriangles].SetIndiceIntersectionPossible(0);
778 CpteurTabTriangles++;
780 PntInit++;//Pass to the next point
782 PntInit++;//Pass the last point of the column
784 TTriangles.SetNbItems(CpteurTabTriangles);
785 const Standard_Integer FinTT = TTriangles.NbItems();
789 //=======================================================================
790 //function : LinkEdges2Triangles
791 //purpose : fill the edge fields in Triangle object for the
792 // two array of triangles.
793 //=======================================================================
794 void IntPolyh_MaillageAffinage::LinkEdges2Triangles()
796 const Standard_Integer FinTT1 = TTriangles1.NbItems();
797 const Standard_Integer FinTT2 = TTriangles2.NbItems();
799 for(Standard_Integer uiui1=0; uiui1<FinTT1; uiui1++) {
800 IntPolyh_Triangle & MyTriangle1=TTriangles1[uiui1];
801 if ( (MyTriangle1.FirstEdge()) == -1 ) {
802 MyTriangle1.SetEdgeandOrientation(1,TEdges1);
803 MyTriangle1.SetEdgeandOrientation(2,TEdges1);
804 MyTriangle1.SetEdgeandOrientation(3,TEdges1);
807 for(Standard_Integer uiui2=0; uiui2<FinTT2; uiui2++) {
808 IntPolyh_Triangle & MyTriangle2=TTriangles2[uiui2];
809 if ( (MyTriangle2.FirstEdge()) == -1 ) {
810 MyTriangle2.SetEdgeandOrientation(1,TEdges2);
811 MyTriangle2.SetEdgeandOrientation(2,TEdges2);
812 MyTriangle2.SetEdgeandOrientation(3,TEdges2);
816 //=======================================================================
817 //function : CommonPartRefinement
818 //purpose : Refine systematicaly all marked triangles of both surfaces
819 // REFINING OF THE COMMON
820 //=======================================================================
821 void IntPolyh_MaillageAffinage::CommonPartRefinement()
823 Standard_Integer FinInit1 = TTriangles1.NbItems();
824 for(Standard_Integer i=0; i<FinInit1; i++) {
825 if(TTriangles1[i].IndiceIntersectionPossible()!=0)
826 TTriangles1[i].MiddleRefinement(i,MaSurface1,TPoints1,TTriangles1,TEdges1);
829 Standard_Integer FinInit2=TTriangles2.NbItems();
830 for(Standard_Integer ii=0; ii<FinInit2; ii++) {
831 if(TTriangles2[ii].IndiceIntersectionPossible()!=0)
832 TTriangles2[ii].MiddleRefinement(ii,MaSurface2,TPoints2,TTriangles2,TEdges2);
836 //=======================================================================
837 //function : LocalSurfaceRefinement
838 //purpose : Refine systematicaly all marked triangles of ONE surface
839 //=======================================================================
840 void IntPolyh_MaillageAffinage::LocalSurfaceRefinement(const Standard_Integer SurfID) {
841 //refine locally, but systematically the chosen surface
843 const Standard_Integer FinInit1 = TTriangles1.NbItems();
844 for(Standard_Integer i=0; i<FinInit1; i++) {
845 if(TTriangles1[i].IndiceIntersectionPossible()!=0)
846 TTriangles1[i].MiddleRefinement(i,MaSurface1,TPoints1,TTriangles1,TEdges1);
851 const Standard_Integer FinInit2 = TTriangles2.NbItems();
852 for(Standard_Integer ii=0; ii<FinInit2; ii++) {
853 if(TTriangles2[ii].IndiceIntersectionPossible()!=0)
854 TTriangles2[ii].MiddleRefinement(ii,MaSurface2,TPoints2,TTriangles2,TEdges2);
858 //=======================================================================
859 //function : ComputeDeflections
860 //purpose : Compute deflection for all triangles of one
861 // surface,and sort min and max of deflections
863 // Calculation of the deflection of all triangles
864 // --> deflection max
865 // --> deflection min
866 //=======================================================================
867 void IntPolyh_MaillageAffinage::ComputeDeflections
868 (const Standard_Integer SurfID)
870 Handle(Adaptor3d_HSurface) MaSurface=(SurfID==1)? MaSurface1:MaSurface2;
871 IntPolyh_ArrayOfPoints &TPoints=(SurfID==1)? TPoints1:TPoints2;
872 IntPolyh_ArrayOfTriangles &TTriangles=(SurfID==1)? TTriangles1:TTriangles2;
873 Standard_Real &FlecheMin=(SurfID==1)? FlecheMin1:FlecheMin2;
874 Standard_Real &FlecheMoy=(SurfID==1)? FlecheMoy1:FlecheMoy2;
875 Standard_Real &FlecheMax=(SurfID==1)? FlecheMax1:FlecheMax2;
877 Standard_Integer CpteurTabFleche=0;
878 FlecheMax=-RealLast();
879 FlecheMin=RealLast();
881 const Standard_Integer FinTT = TTriangles.NbItems();
883 for(CpteurTabFleche=0; CpteurTabFleche<FinTT; CpteurTabFleche++) {
884 IntPolyh_Triangle &Triangle = TTriangles[CpteurTabFleche];
885 if ( Triangle.GetFleche() < 0) { //pas normal
889 Triangle.TriangleDeflection(MaSurface, TPoints);
890 Standard_Real Fleche=Triangle.GetFleche();
892 if (Fleche > FlecheMax)
894 if (Fleche < FlecheMin)
899 //=======================================================================
900 //function : TrianglesDeflectionsRefinementBSB
901 //purpose : Refine both surfaces using BoundSortBox as --
902 // rejection. The criterions used to refine a --
903 // triangle are: The deflection The size of the --
904 // bounding boxes (one surface may be very small
905 // compared to the other)
906 //=======================================================================
907 void IntPolyh_MaillageAffinage::TrianglesDeflectionsRefinementBSB()
909 const Standard_Integer FinTT1 = TTriangles1.NbItems();
910 const Standard_Integer FinTT2 = TTriangles2.NbItems();
912 ComputeDeflections(1);
913 // To estimate a surface in general it can be interesting
914 //to calculate all deflections
915 //-- Check deflection at output
917 Standard_Real FlecheCritique1;
918 if(FlecheMin1>FlecheMax1) {
921 else {//fleche min + (flechemax-flechemin) * 80/100
922 FlecheCritique1 = FlecheMin1*0.2+FlecheMax1*0.8;
925 ComputeDeflections(2);
926 //-- Check arrows at output
928 Standard_Real FlecheCritique2;
929 if(FlecheMin2>FlecheMax2) {
933 else {//fleche min + (flechemax-flechemin) * 80/100
934 FlecheCritique2 = FlecheMin2*0.2+FlecheMax2*0.8;
938 Bnd_BoundSortBox BndBSB;
939 Standard_Real diag1,diag2;
940 Standard_Real x0,y0,z0,x1,y1,z1;
942 //The greatest of two bounding boxes created in FillArrayOfPoints is found.
943 //Then this value is weighted depending on the discretization
944 //(NbSamplesU and NbSamplesV)
945 MyBox1.Get(x0,y0,z0,x1,y1,z1);
946 x0-=x1; y0-=y1; z0-=z1;
947 diag1=x0*x0+y0*y0+z0*z0;
948 const Standard_Real NbSamplesUV1=Standard_Real(NbSamplesU1) * Standard_Real(NbSamplesV1);
951 MyBox2.Get(x0,y0,z0,x1,y1,z1);
952 x0-=x1; y0-=y1; z0-=z1;
953 diag2=x0*x0+y0*y0+z0*z0;
954 const Standard_Real NbSamplesUV2=Standard_Real(NbSamplesU2) * Standard_Real(NbSamplesV2);
957 //-- The surface with the greatest bounding box is "discretized"
959 //Standard_Integer NbInterTentees=0;
963 if(FlecheCritique2<diag1) {//the corresponding sizes are not too disproportional
965 Handle(Bnd_HArray1OfBox) HBnd = new Bnd_HArray1OfBox(1,FinTT2);
967 for(Standard_Integer i=0; i<FinTT2; i++){
968 if (TTriangles2[i].IndiceIntersectionPossible()!=0) {
970 const IntPolyh_Triangle& T=TTriangles2[i];
971 const IntPolyh_Point& PA=TPoints2[T.FirstPoint()];
972 const IntPolyh_Point& PB=TPoints2[T.SecondPoint()];
973 const IntPolyh_Point& PC=TPoints2[T.ThirdPoint()];
974 gp_Pnt pntA(PA.X(),PA.Y(),PA.Z());
975 gp_Pnt pntB(PB.X(),PB.Y(),PB.Z());
976 gp_Pnt pntC(PC.X(),PC.Y(),PC.Z());
977 b.Add(pntA);//Box b, which contains triangle i of surface 2 is created./
980 b.Enlarge(T.GetFleche()+MyTolerance);
981 HBnd->SetValue(i+1,b);//Box b is added in the array HBnd
985 //Inititalization of the boundary, sorting of boxes
986 BndBSB.Initialize(HBnd);//contains boxes of 2
988 Standard_Integer FinTT1Init=FinTT1;
989 for(Standard_Integer i_S1=0; i_S1<FinTT1Init; i_S1++) {
990 if(TTriangles1[i_S1].IndiceIntersectionPossible()!=0) {
991 //-- Loop on the boxes of mesh 1
993 const IntPolyh_Triangle& T=TTriangles1[i_S1];
994 const IntPolyh_Point& PA=TPoints1[T.FirstPoint()];
995 const IntPolyh_Point& PB=TPoints1[T.SecondPoint()];
996 const IntPolyh_Point& PC=TPoints1[T.ThirdPoint()];
997 gp_Pnt pntA(PA.X(),PA.Y(),PA.Z());
998 gp_Pnt pntB(PB.X(),PB.Y(),PB.Z());
999 gp_Pnt pntC(PC.X(),PC.Y(),PC.Z());
1003 b.Enlarge(T.GetFleche());
1004 //-- List of boxes of 2, which touch this box (of 1)
1005 const TColStd_ListOfInteger& ListeOf2 = BndBSB.Compare(b);
1007 if((ListeOf2.IsEmpty())==0) {
1008 IntPolyh_Triangle &Triangle1 = TTriangles1[i_S1];
1009 if(Triangle1.GetFleche()>FlecheCritique1)
1010 Triangle1.MiddleRefinement(i_S1,MaSurface1,TPoints1,
1011 TTriangles1, TEdges1);
1013 for (TColStd_ListIteratorOfListOfInteger Iter(ListeOf2);
1016 Standard_Integer i_S2=Iter.Value()-1;
1017 //if the box of s1 contacts with the boxes of s2
1018 //the arrow of the triangle is checked
1019 IntPolyh_Triangle & Triangle2 = TTriangles2[i_S2];
1020 if(Triangle2.IndiceIntersectionPossible()!=0)
1021 if(Triangle2.GetFleche()>FlecheCritique2)
1022 Triangle2.MiddleRefinement( i_S2, MaSurface2, TPoints2,
1023 TTriangles2, TEdges2);
1030 //--------------------------------------------------------------------
1031 //FlecheCritique2 > diag1
1035 Handle(Bnd_HArray1OfBox) HBnd = new Bnd_HArray1OfBox(1,FinTT2);
1037 for(Standard_Integer i=0; i<FinTT2; i++){
1038 if (TTriangles2[i].IndiceIntersectionPossible()!=0) {
1040 const IntPolyh_Triangle& T=TTriangles2[i];
1041 const IntPolyh_Point& PA=TPoints2[T.FirstPoint()];
1042 const IntPolyh_Point& PB=TPoints2[T.SecondPoint()];
1043 const IntPolyh_Point& PC=TPoints2[T.ThirdPoint()];
1044 gp_Pnt pntA(PA.X(),PA.Y(),PA.Z());
1045 gp_Pnt pntB(PB.X(),PB.Y(),PB.Z());
1046 gp_Pnt pntC(PC.X(),PC.Y(),PC.Z());
1047 b.Add(pntA);//Box b, which contains triangle i of surface 2 is created/
1050 b.Enlarge(T.GetFleche()+MyTolerance);
1051 //-- BndBSB.Add(b,i+1);
1052 HBnd->SetValue(i+1,b);//Box b is added in array HBnd
1056 //Inititalization of the ouput bounding box
1057 BndBSB.Initialize(HBnd);//contains boxes of 2
1060 //The bounding box Be1 of surface1 is compared BSB of surface2
1061 const TColStd_ListOfInteger& ListeOf2 = BndBSB.Compare(MyBox1);
1063 if((ListeOf2.IsEmpty())==0) {
1064 //if the bounding box Be1 of s1 contacts with
1065 //the boxes of s2 the deflection of triangle of s2 is checked
1067 // Be1 is very small in relation to Be2
1068 //The criterion of refining for surface2 depends on the size of Be1
1069 //As it is known that this criterion should be minimized,
1070 //the smallest side of the bounding box is taken
1071 Standard_Real x0,x1,y0,y1,z0,z1;
1072 MyBox1.Get(x0,y0,z0,x1,y1,z1);
1073 Standard_Real dx=Abs(x1-x0);
1074 Standard_Real dy=Abs(y1-y0);
1075 Standard_Real diag=Abs(z1-z0);
1076 Standard_Real dd=-1.0;
1081 if (diag>dd) diag=dd;
1083 //if Be1 contacts with the boxes of s2, the deflection
1084 //of the triangles of s2 is checked (greater)
1085 //in relation to the size of Be1 (smaller)
1086 for (TColStd_ListIteratorOfListOfInteger Iter(ListeOf2);
1089 Standard_Integer i_S2=Iter.Value()-1;
1091 IntPolyh_Triangle & Triangle2=TTriangles2[i_S2];
1092 if(Triangle2.IndiceIntersectionPossible()) {
1094 //calculation of the criterion of refining
1095 //The deflection of the greater is compared to the size of the smaller
1096 Standard_Real CritereAffinage=0.0;
1097 Standard_Real DiagPonderation=0.5;
1098 CritereAffinage = diag*DiagPonderation;
1099 if(Triangle2.GetFleche()>CritereAffinage)
1100 Triangle2.MultipleMiddleRefinement2(CritereAffinage, MyBox1, i_S2,
1101 MaSurface2, TPoints2,
1102 TTriangles2,TEdges2);
1104 else Triangle2.MiddleRefinement(i_S2,MaSurface2,TPoints2,
1105 TTriangles2, TEdges2);
1113 else { //-- The greater is discretised
1115 if(FlecheCritique1<diag2) {//the respective sizes are not to much disproportional
1117 Handle(Bnd_HArray1OfBox) HBnd = new Bnd_HArray1OfBox(1,FinTT1);
1119 for(Standard_Integer i=0; i<FinTT1; i++){
1120 if(TTriangles1[i].IndiceIntersectionPossible()!=0) {
1122 const IntPolyh_Triangle& T=TTriangles1[i];
1123 const IntPolyh_Point& PA=TPoints1[T.FirstPoint()];
1124 const IntPolyh_Point& PB=TPoints1[T.SecondPoint()];
1125 const IntPolyh_Point& PC=TPoints1[T.ThirdPoint()];
1126 gp_Pnt pntA(PA.X(),PA.Y(),PA.Z());
1127 gp_Pnt pntB(PB.X(),PB.Y(),PB.Z());
1128 gp_Pnt pntC(PC.X(),PC.Y(),PC.Z());
1129 b.Add(pntA);//Box b, which contains triangle i of surface 2 is created.
1132 b.Enlarge(T.GetFleche()+MyTolerance);
1133 HBnd->SetValue(i+1,b);//Boite b is added in the array HBnd
1136 BndBSB.Initialize(HBnd);
1138 Standard_Integer FinTT2init=FinTT2;
1139 for(Standard_Integer i_S2=0; i_S2<FinTT2init; i_S2++) {
1140 if (TTriangles2[i_S2].IndiceIntersectionPossible()!=0) {
1141 //-- Loop on the boxes of mesh 2
1143 const IntPolyh_Triangle& T=TTriangles2[i_S2];
1144 const IntPolyh_Point& PA=TPoints2[T.FirstPoint()];
1145 const IntPolyh_Point& PB=TPoints2[T.SecondPoint()];
1146 const IntPolyh_Point& PC=TPoints2[T.ThirdPoint()];
1147 gp_Pnt pntA(PA.X(),PA.Y(),PA.Z());
1148 gp_Pnt pntB(PB.X(),PB.Y(),PB.Z());
1149 gp_Pnt pntC(PC.X(),PC.Y(),PC.Z());
1153 b.Enlarge(T.GetFleche()+MyTolerance);
1154 //-- List of boxes of 1 touching this box (of 2)
1155 const TColStd_ListOfInteger& ListeOf1 = BndBSB.Compare(b);
1156 IntPolyh_Triangle & Triangle2=TTriangles2[i_S2];
1157 if((ListeOf1.IsEmpty())==0) {
1159 if(Triangle2.GetFleche()>FlecheCritique2)
1160 Triangle2.MiddleRefinement(i_S2,MaSurface2,TPoints2,
1161 TTriangles2, TEdges2);
1163 for (TColStd_ListIteratorOfListOfInteger Iter(ListeOf1);
1166 Standard_Integer i_S1=Iter.Value()-1;
1167 IntPolyh_Triangle & Triangle1=TTriangles1[i_S1];
1168 if (Triangle1.IndiceIntersectionPossible())
1169 if(Triangle1.GetFleche()>FlecheCritique1)
1170 Triangle1.MiddleRefinement(i_S1,MaSurface1,TPoints1,
1171 TTriangles1, TEdges1);
1177 //-----------------------------------------------------------------------------
1178 else {// FlecheCritique1>diag2
1181 Handle(Bnd_HArray1OfBox) HBnd = new Bnd_HArray1OfBox(1,FinTT1);
1183 for(Standard_Integer i=0; i<FinTT1; i++){
1184 if (TTriangles1[i].IndiceIntersectionPossible()!=0) {
1186 const IntPolyh_Triangle& T=TTriangles1[i];
1187 const IntPolyh_Point& PA=TPoints1[T.FirstPoint()];
1188 const IntPolyh_Point& PB=TPoints1[T.SecondPoint()];
1189 const IntPolyh_Point& PC=TPoints1[T.ThirdPoint()];
1190 gp_Pnt pntA(PA.X(),PA.Y(),PA.Z());
1191 gp_Pnt pntB(PB.X(),PB.Y(),PB.Z());
1192 gp_Pnt pntC(PC.X(),PC.Y(),PC.Z());
1193 b.Add(pntA);//Box b, which contains triangle i of surface 1 is created./
1196 b.Enlarge(T.GetFleche()+MyTolerance);
1197 HBnd->SetValue(i+1,b);//Box b is added in the array HBnd
1201 //Inititalisation of the boundary output box
1202 BndBSB.Initialize(HBnd);//contains boxes of 1
1204 //Bounding box Be2 of surface2 is compared to BSB of surface1
1205 const TColStd_ListOfInteger& ListeOf1 = BndBSB.Compare(MyBox2);
1207 if((ListeOf1.IsEmpty())==0) {
1208 //if the bounding box Be2 of s2 contacts
1209 //with boxes of s1 the deflection of the triangle of s1 is checked
1211 // Be2 is very small compared to Be1
1212 //The criterion of refining for surface1 depends on the size of Be2
1213 //As this criterion should be minimized,
1214 //the smallest side of the bounding box is taken
1215 Standard_Real x0,x1,y0,y1,z0,z1;
1216 MyBox2.Get(x0,y0,z0,x1,y1,z1);
1217 Standard_Real dx=Abs(x1-x0);
1218 Standard_Real dy=Abs(y1-y0);
1219 Standard_Real diag=Abs(z1-z0);
1220 Standard_Real dd=-1.0;
1225 if (diag>dd) diag=dd;
1227 //if Be2 contacts with boxes of s1, the deflection of
1228 //triangles of s1 (greater) is checked
1229 //comparatively to the size of Be2 (smaller).
1230 for (TColStd_ListIteratorOfListOfInteger Iter(ListeOf1);
1233 Standard_Integer i_S1=Iter.Value()-1;
1235 IntPolyh_Triangle & Triangle1=TTriangles1[i_S1];
1236 if(Triangle1.IndiceIntersectionPossible()) {
1238 //calculation of the criterion of refining
1239 //The deflection of the greater is compared
1240 //with the size of the smaller.
1241 Standard_Real CritereAffinage=0.0;
1242 Standard_Real DiagPonderation=0.5;
1243 CritereAffinage = diag*DiagPonderation;;
1244 if(Triangle1.GetFleche()>CritereAffinage)
1245 Triangle1.MultipleMiddleRefinement2(CritereAffinage,MyBox2, i_S1,
1246 MaSurface1, TPoints1,
1247 TTriangles1, TEdges1);
1249 else Triangle1.MiddleRefinement(i_S1,MaSurface1,TPoints1,
1250 TTriangles1, TEdges1);
1258 //=======================================================================
1260 //purpose : This function is used for the function project6
1261 //=======================================================================
1262 inline Standard_Real maxSR(const Standard_Real a,
1263 const Standard_Real b,
1264 const Standard_Real c)
1266 Standard_Real t = a;
1271 //=======================================================================
1273 //purpose : This function is used for the function project6
1274 //=======================================================================
1275 inline Standard_Real minSR(const Standard_Real a,
1276 const Standard_Real b,
1277 const Standard_Real c)
1279 Standard_Real t = a;
1284 //=======================================================================
1285 //function : project6
1286 //purpose : This function is used for the function TriContact
1287 //=======================================================================
1288 Standard_Integer project6(const IntPolyh_Point &ax,
1289 const IntPolyh_Point &p1,
1290 const IntPolyh_Point &p2,
1291 const IntPolyh_Point &p3,
1292 const IntPolyh_Point &q1,
1293 const IntPolyh_Point &q2,
1294 const IntPolyh_Point &q3)
1296 Standard_Real P1 = ax.Dot(p1);
1297 Standard_Real P2 = ax.Dot(p2);
1298 Standard_Real P3 = ax.Dot(p3);
1299 Standard_Real Q1 = ax.Dot(q1);
1300 Standard_Real Q2 = ax.Dot(q2);
1301 Standard_Real Q3 = ax.Dot(q3);
1303 Standard_Real mx1 = maxSR(P1, P2, P3);
1304 Standard_Real mn1 = minSR(P1, P2, P3);
1305 Standard_Real mx2 = maxSR(Q1, Q2, Q3);
1306 Standard_Real mn2 = minSR(Q1, Q2, Q3);
1308 if (mn1 > mx2) return 0;
1309 if (mn2 > mx1) return 0;
1312 //=======================================================================
1313 //function : TriContact
1314 //purpose : This fonction Check if two triangles are in
1315 // contact or no, return 1 if yes, return 0
1317 //=======================================================================
1318 Standard_Integer IntPolyh_MaillageAffinage::TriContact
1319 (const IntPolyh_Point &P1,
1320 const IntPolyh_Point &P2,
1321 const IntPolyh_Point &P3,
1322 const IntPolyh_Point &Q1,
1323 const IntPolyh_Point &Q2,
1324 const IntPolyh_Point &Q3,
1325 Standard_Real &Angle) const
1328 The first triangle is (p1,p2,p3). The other is (q1,q2,q3).
1329 The edges are (e1,e2,e3) and (f1,f2,f3).
1330 The normals are n1 and m1
1331 Outwards are (g1,g2,g3) and (h1,h2,h3).*/
1333 IntPolyh_Point p1, p2, p3;
1334 IntPolyh_Point q1, q2, q3;
1335 IntPolyh_Point e1, e2, e3;
1336 IntPolyh_Point f1, f2, f3;
1337 IntPolyh_Point g1, g2, g3;
1338 IntPolyh_Point h1, h2, h3;
1339 IntPolyh_Point n1, m1;
1342 IntPolyh_Point ef11, ef12, ef13;
1343 IntPolyh_Point ef21, ef22, ef23;
1344 IntPolyh_Point ef31, ef32, ef33;
1346 z.SetX(0.0); z.SetY(0.0); z.SetZ(0.0);
1348 if(maxSR(P1.X(),P2.X(),P3.X())<minSR(Q1.X(),Q2.X(),Q3.X())) return(0);
1349 if(maxSR(P1.Y(),P2.Y(),P3.Y())<minSR(Q1.Y(),Q2.Y(),Q3.Y())) return(0);
1350 if(maxSR(P1.Z(),P2.Z(),P3.Z())<minSR(Q1.Z(),Q2.Z(),Q3.Z())) return(0);
1352 if(minSR(P1.X(),P2.X(),P3.X())>maxSR(Q1.X(),Q2.X(),Q3.X())) return(0);
1353 if(minSR(P1.Y(),P2.Y(),P3.Y())>maxSR(Q1.Y(),Q2.Y(),Q3.Y())) return(0);
1354 if(minSR(P1.Z(),P2.Z(),P3.Z())>maxSR(Q1.Z(),Q2.Z(),Q3.Z())) return(0);
1356 p1.SetX(P1.X() - P1.X()); p1.SetY(P1.Y() - P1.Y()); p1.SetZ(P1.Z() - P1.Z());
1357 p2.SetX(P2.X() - P1.X()); p2.SetY(P2.Y() - P1.Y()); p2.SetZ(P2.Z() - P1.Z());
1358 p3.SetX(P3.X() - P1.X()); p3.SetY(P3.Y() - P1.Y()); p3.SetZ(P3.Z() - P1.Z());
1360 q1.SetX(Q1.X() - P1.X()); q1.SetY(Q1.Y() - P1.Y()); q1.SetZ(Q1.Z() - P1.Z());
1361 q2.SetX(Q2.X() - P1.X()); q2.SetY(Q2.Y() - P1.Y()); q2.SetZ(Q2.Z() - P1.Z());
1362 q3.SetX(Q3.X() - P1.X()); q3.SetY(Q3.Y() - P1.Y()); q3.SetZ(Q3.Z() - P1.Z());
1364 e1.SetX(p2.X() - p1.X()); e1.SetY(p2.Y() - p1.Y()); e1.SetZ(p2.Z() - p1.Z());
1365 e2.SetX(p3.X() - p2.X()); e2.SetY(p3.Y() - p2.Y()); e2.SetZ(p3.Z() - p2.Z());
1366 e3.SetX(p1.X() - p3.X()); e3.SetY(p1.Y() - p3.Y()); e3.SetZ(p1.Z() - p3.Z());
1368 f1.SetX(q2.X() - q1.X()); f1.SetY(q2.Y() - q1.Y()); f1.SetZ(q2.Z() - q1.Z());
1369 f2.SetX(q3.X() - q2.X()); f2.SetY(q3.Y() - q2.Y()); f2.SetZ(q3.Z() - q2.Z());
1370 f3.SetX(q1.X() - q3.X()); f3.SetY(q1.Y() - q3.Y()); f3.SetZ(q1.Z() - q3.Z());
1372 n1.Cross(e1, e2); //normal to the first triangle
1373 m1.Cross(f1, f2); //normal to the second triangle
1392 // Now the testing is done
1394 if (!project6(n1, p1, p2, p3, q1, q2, q3)) return 0; //T2 is not higher or lower than T1
1395 if (!project6(m1, p1, p2, p3, q1, q2, q3)) return 0; //T1 is not higher of lower than T2
1397 if (!project6(ef11, p1, p2, p3, q1, q2, q3)) return 0;
1398 if (!project6(ef12, p1, p2, p3, q1, q2, q3)) return 0;
1399 if (!project6(ef13, p1, p2, p3, q1, q2, q3)) return 0;
1400 if (!project6(ef21, p1, p2, p3, q1, q2, q3)) return 0;
1401 if (!project6(ef22, p1, p2, p3, q1, q2, q3)) return 0;
1402 if (!project6(ef23, p1, p2, p3, q1, q2, q3)) return 0;
1403 if (!project6(ef31, p1, p2, p3, q1, q2, q3)) return 0;
1404 if (!project6(ef32, p1, p2, p3, q1, q2, q3)) return 0;
1405 if (!project6(ef33, p1, p2, p3, q1, q2, q3)) return 0;
1407 if (!project6(g1, p1, p2, p3, q1, q2, q3)) return 0; //T2 is outside of T1 in the plane of T1
1408 if (!project6(g2, p1, p2, p3, q1, q2, q3)) return 0; //T2 is outside of T1 in the plane of T1
1409 if (!project6(g3, p1, p2, p3, q1, q2, q3)) return 0; //T2 is outside of T1 in the plane of T1
1410 if (!project6(h1, p1, p2, p3, q1, q2, q3)) return 0; //T1 is outside of T2 in the plane of T2
1411 if (!project6(h2, p1, p2, p3, q1, q2, q3)) return 0; //T1 is outside of T2 in the plane of T2
1412 if (!project6(h3, p1, p2, p3, q1, q2, q3)) return 0; //T1 is outside of T2 in the plane of T2
1414 //Calculation of cosinus angle between two normals
1415 Standard_Real SqModn1=-1.0;
1416 Standard_Real SqModm1=-1.0;
1417 SqModn1=n1.SquareModulus();
1418 if (SqModn1>SquareMyConfusionPrecision){
1419 SqModm1=m1.SquareModulus();
1421 if (SqModm1>SquareMyConfusionPrecision) {
1422 Angle=(n1.Dot(m1))/(sqrt(SqModn1)*sqrt(SqModm1));
1426 //=======================================================================
1427 //function : TestNbPoints
1428 //purpose : This function is used by StartingPointsResearch() to control
1429 // the number of points found keep the result in conformity (1 or 2 points)
1430 // void TestNbPoints(const Standard_Integer TriSurfID,
1431 //=======================================================================
1432 void TestNbPoints(const Standard_Integer ,
1433 Standard_Integer &NbPoints,
1434 Standard_Integer &NbPointsTotal,
1435 const IntPolyh_StartPoint &Pt1,
1436 const IntPolyh_StartPoint &Pt2,
1437 IntPolyh_StartPoint &SP1,
1438 IntPolyh_StartPoint &SP2)
1440 // already checked in TriangleEdgeContact2
1441 // if( (NbPoints==2)&&(Pt1.CheckSameSP(Pt2)) ) NbPoints=1;
1447 if ( (NbPoints==1)&&(NbPointsTotal==0) ) {
1451 else if ( (NbPoints==1)&&(NbPointsTotal==1) ) {
1452 if(Pt1.CheckSameSP(SP1)!=1) {
1457 else if( (NbPoints==1)&&(NbPointsTotal==2) ) {
1458 if ( (SP1.CheckSameSP(Pt1))||(SP2.CheckSameSP(Pt1)) )
1460 else NbPointsTotal=3;
1462 else if( (NbPoints==2)&&(NbPointsTotal==0) ) {
1467 else if( (NbPoints==2)&&(NbPointsTotal==1) ) {//there is also Pt1 != Pt2
1468 if(SP1.CheckSameSP(Pt1)) {
1472 else if (SP1.CheckSameSP(Pt2)) {
1476 else NbPointsTotal=3;///case SP1!=Pt1 && SP1!=Pt2!
1478 else if( (NbPoints==2)&&(NbPointsTotal==2) ) {//there is also SP1!=SP2
1479 if( (SP1.CheckSameSP(Pt1))||(SP1.CheckSameSP(Pt2)) ) {
1480 if( (SP2.CheckSameSP(Pt1))||(SP2.CheckSameSP(Pt2)) )
1482 else NbPointsTotal=3;
1484 else NbPointsTotal=3;
1488 //=======================================================================
1489 //function : StartingPointsResearch
1491 //=======================================================================
1492 Standard_Integer IntPolyh_MaillageAffinage::StartingPointsResearch
1493 (const Standard_Integer T1,
1494 const Standard_Integer T2,
1495 IntPolyh_StartPoint &SP1,
1496 IntPolyh_StartPoint &SP2) const
1498 const IntPolyh_Point &P1=TPoints1[TTriangles1[T1].FirstPoint()];
1499 const IntPolyh_Point &P2=TPoints1[TTriangles1[T1].SecondPoint()];
1500 const IntPolyh_Point &P3=TPoints1[TTriangles1[T1].ThirdPoint()];
1501 const IntPolyh_Point &Q1=TPoints2[TTriangles2[T2].FirstPoint()];
1502 const IntPolyh_Point &Q2=TPoints2[TTriangles2[T2].SecondPoint()];
1503 const IntPolyh_Point &Q3=TPoints2[TTriangles2[T2].ThirdPoint()];
1506 /* The first triangle is (p1,p2,p3). The other is (q1,q2,q3).
1507 The sides are (e1,e2,e3) and (f1,f2,f3).
1508 The normals are n1 and m1*/
1510 const IntPolyh_Point e1=P2-P1;
1511 const IntPolyh_Point e2=P3-P2;
1512 const IntPolyh_Point e3=P1-P3;
1514 const IntPolyh_Point f1=Q2-Q1;
1515 const IntPolyh_Point f2=Q3-Q2;
1516 const IntPolyh_Point f3=Q1-Q3;
1519 IntPolyh_Point nn1,mm1;
1520 nn1.Cross(e1, e2); //normal of the first triangle
1521 mm1.Cross(f1, f2); //normal of the second triangle
1523 Standard_Real nn1modulus, mm1modulus;
1524 nn1modulus=sqrt(nn1.SquareModulus());
1525 mm1modulus=sqrt(mm1.SquareModulus());
1527 //-------------------------------------------------------
1528 ///calculation of intersection points between two triangles
1529 //-------------------------------------------------------
1530 Standard_Integer NbPoints=0;
1531 Standard_Integer NbPointsTotal=0;
1532 IntPolyh_StartPoint Pt1,Pt2;
1536 if(Abs(nn1modulus)<MyConfusionPrecision){//10.0e-20) {
1540 const IntPolyh_Point n1=nn1.Divide(nn1modulus);
1542 if(NbPointsTotal<2) {
1543 NbPoints=TriangleEdgeContact(1,1,P1,P2,P3,e1,e2,e3,Q1,Q2,f1,n1,Pt1,Pt2);
1544 TestNbPoints(1,NbPoints,NbPointsTotal,Pt1,Pt2,SP1,SP2);
1547 if(NbPointsTotal<2) {
1548 NbPoints=TriangleEdgeContact(1,2,P1,P2,P3,e1,e2,e3,Q2,Q3,f2,n1,Pt1,Pt2);
1549 TestNbPoints(1,NbPoints,NbPointsTotal,Pt1,Pt2,SP1,SP2);
1552 if(NbPointsTotal<2) {
1553 NbPoints=TriangleEdgeContact(1,3,P1,P2,P3,e1,e2,e3,Q3,Q1,f3,n1,Pt1,Pt2);
1554 TestNbPoints(1,NbPoints,NbPointsTotal,Pt1,Pt2,SP1,SP2);
1559 if(Abs(mm1modulus)<MyConfusionPrecision) {//10.0e-20){
1563 const IntPolyh_Point m1=mm1.Divide(mm1modulus);
1565 if(NbPointsTotal<2) {
1566 NbPoints=TriangleEdgeContact(2,1,Q1,Q2,Q3,f1,f2,f3,P1,P2,e1,m1,Pt1,Pt2);
1567 TestNbPoints(2,NbPoints,NbPointsTotal,Pt1,Pt2,SP1,SP2);
1570 if(NbPointsTotal<2) {
1571 NbPoints=TriangleEdgeContact(2,2,Q1,Q2,Q3,f1,f2,f3,P2,P3,e2,m1,Pt1,Pt2);
1572 TestNbPoints(2,NbPoints,NbPointsTotal,Pt1,Pt2,SP1,SP2);
1575 if(NbPointsTotal<2) {
1576 NbPoints=TriangleEdgeContact(2,3,Q1,Q2,Q3,f1,f2,f3,P3,P1,e3,m1,Pt1,Pt2);
1577 TestNbPoints(2,NbPoints,NbPointsTotal,Pt1,Pt2,SP1,SP2);
1581 /* if( (NbPointsTotal >1)&&( Abs(SP1.U1()-SP2.U1())<MyConfusionPrecision)
1582 &&(Abs(SP1.V1()-SP2.V1())<MyConfusionPrecision) )*/
1583 if( (NbPoints)&&(SP1.CheckSameSP(SP2)) )
1586 SP1.SetCoupleValue(T1,T2);
1587 SP2.SetCoupleValue(T1,T2);
1588 return (NbPointsTotal);
1590 //=======================================================================
1591 //function : StartingPointsResearch2
1592 //purpose : From two triangles compute intersection points.
1593 // If I found more than two intersection points
1594 // it means that those triangle are coplanar
1595 //=======================================================================
1596 Standard_Integer IntPolyh_MaillageAffinage::StartingPointsResearch2
1597 (const Standard_Integer T1,
1598 const Standard_Integer T2,
1599 IntPolyh_StartPoint &SP1,
1600 IntPolyh_StartPoint &SP2) const
1602 const IntPolyh_Triangle &Tri1=TTriangles1[T1];
1603 const IntPolyh_Triangle &Tri2=TTriangles2[T2];
1605 const IntPolyh_Point &P1=TPoints1[Tri1.FirstPoint()];
1606 const IntPolyh_Point &P2=TPoints1[Tri1.SecondPoint()];
1607 const IntPolyh_Point &P3=TPoints1[Tri1.ThirdPoint()];
1608 const IntPolyh_Point &Q1=TPoints2[Tri2.FirstPoint()];
1609 const IntPolyh_Point &Q2=TPoints2[Tri2.SecondPoint()];
1610 const IntPolyh_Point &Q3=TPoints2[Tri2.ThirdPoint()];
1614 /* The first triangle is (p1,p2,p3). The other is (q1,q2,q3).
1615 The sides are (e1,e2,e3) and (f1,f2,f3).
1616 The normals are n1 and m1*/
1618 const IntPolyh_Point e1=P2-P1;
1619 const IntPolyh_Point e2=P3-P2;
1620 const IntPolyh_Point e3=P1-P3;
1622 const IntPolyh_Point f1=Q2-Q1;
1623 const IntPolyh_Point f2=Q3-Q2;
1624 const IntPolyh_Point f3=Q1-Q3;
1627 IntPolyh_Point nn1,mm1;
1628 nn1.Cross(e1, e2); //normal to the first triangle
1629 mm1.Cross(f1, f2); //normal to the second triangle
1631 Standard_Real nn1modulus, mm1modulus;
1632 nn1modulus=sqrt(nn1.SquareModulus());
1633 mm1modulus=sqrt(mm1.SquareModulus());
1635 //-------------------------------------------------
1636 ///calculation of intersections points between triangles
1637 //-------------------------------------------------
1638 Standard_Integer NbPoints=0;
1639 Standard_Integer NbPointsTotal=0;
1643 if(Abs(nn1modulus)<MyConfusionPrecision) {//10.0e-20){
1647 const IntPolyh_Point n1=nn1.Divide(nn1modulus);
1649 if(NbPointsTotal<3) {
1650 IntPolyh_StartPoint Pt1,Pt2;
1651 NbPoints=TriangleEdgeContact2(1,1,Tri1,Tri2,P1,P2,P3,e1,e2,e3,Q1,Q2,f1,n1,Pt1,Pt2);
1652 TestNbPoints(1,NbPoints,NbPointsTotal,Pt1,Pt2,SP1,SP2);
1655 if(NbPointsTotal<3) {
1656 IntPolyh_StartPoint Pt1,Pt2;
1657 NbPoints=TriangleEdgeContact2(1,2,Tri1,Tri2,P1,P2,P3,e1,e2,e3,Q2,Q3,f2,n1,Pt1,Pt2);
1658 TestNbPoints(1,NbPoints,NbPointsTotal,Pt1,Pt2,SP1,SP2);
1661 if(NbPointsTotal<3) {
1662 IntPolyh_StartPoint Pt1,Pt2;
1663 NbPoints=TriangleEdgeContact2(1,3,Tri1,Tri2,P1,P2,P3,e1,e2,e3,Q3,Q1,f3,n1,Pt1,Pt2);
1664 TestNbPoints(1,NbPoints,NbPointsTotal,Pt1,Pt2,SP1,SP2);
1669 if(Abs(mm1modulus)<MyConfusionPrecision) {//10.0e-20){
1673 const IntPolyh_Point m1=mm1.Divide(mm1modulus);
1675 if(NbPointsTotal<3) {
1676 IntPolyh_StartPoint Pt1,Pt2;
1677 NbPoints=TriangleEdgeContact2(2,1,Tri1,Tri2,Q1,Q2,Q3,f1,f2,f3,P1,P2,e1,m1,Pt1,Pt2);
1678 TestNbPoints(2,NbPoints,NbPointsTotal,Pt1,Pt2,SP1,SP2);
1681 if(NbPointsTotal<3) {
1682 IntPolyh_StartPoint Pt1,Pt2;
1683 NbPoints=TriangleEdgeContact2(2,2,Tri1,Tri2,Q1,Q2,Q3,f1,f2,f3,P2,P3,e2,m1,Pt1,Pt2);
1684 TestNbPoints(2,NbPoints,NbPointsTotal,Pt1,Pt2,SP1,SP2);
1687 if(NbPointsTotal<3) {
1688 IntPolyh_StartPoint Pt1,Pt2;
1689 NbPoints=TriangleEdgeContact2(2,3,Tri1,Tri2,Q1,Q2,Q3,f1,f2,f3,P3,P1,e3,m1,Pt1,Pt2);
1690 TestNbPoints(2,NbPoints,NbPointsTotal,Pt1,Pt2,SP1,SP2);
1693 // no need already checked in TestNbPoints
1694 // if( (NbPointsTotal==2)&&(SP1.CheckSameSP(SP2)) ) {
1696 //SP1.SetCoupleValue(T1,T2);
1699 if(NbPointsTotal==2) {
1700 SP1.SetCoupleValue(T1,T2);
1701 SP2.SetCoupleValue(T1,T2);
1703 else if(NbPointsTotal==1)
1704 SP1.SetCoupleValue(T1,T2);
1705 else if(NbPointsTotal==3)
1706 SP1.SetCoupleValue(T1,T2);
1708 return (NbPointsTotal);
1710 //=======================================================================
1711 //function : NextStartingPointsResearch
1713 //=======================================================================
1714 Standard_Integer IntPolyh_MaillageAffinage::NextStartingPointsResearch
1715 (const Standard_Integer T1,
1716 const Standard_Integer T2,
1717 const IntPolyh_StartPoint &SPInit,
1718 IntPolyh_StartPoint &SPNext) const
1720 Standard_Integer NbPointsTotal=0;
1721 if( (T1<0)||(T2<0) ) NbPointsTotal=0;
1723 const IntPolyh_Point &P1=TPoints1[TTriangles1[T1].FirstPoint()];
1724 const IntPolyh_Point &P2=TPoints1[TTriangles1[T1].SecondPoint()];
1725 const IntPolyh_Point &P3=TPoints1[TTriangles1[T1].ThirdPoint()];
1726 const IntPolyh_Point &Q1=TPoints2[TTriangles2[T2].FirstPoint()];
1727 const IntPolyh_Point &Q2=TPoints2[TTriangles2[T2].SecondPoint()];
1728 const IntPolyh_Point &Q3=TPoints2[TTriangles2[T2].ThirdPoint()];
1730 /* The first triangle is (p1,p2,p3). The other is (q1,q2,q3).
1731 The sides are (e1,e2,e3) and (f1,f2,f3).
1732 The normals are n1 and m1*/
1734 const IntPolyh_Point e1=P2-P1;
1735 const IntPolyh_Point e2=P3-P2;
1736 const IntPolyh_Point e3=P1-P3;
1738 const IntPolyh_Point f1=Q2-Q1;
1739 const IntPolyh_Point f2=Q3-Q2;
1740 const IntPolyh_Point f3=Q1-Q3;
1742 IntPolyh_Point nn1,mm1;
1743 nn1.Cross(e1, e2); //normal to the first triangle
1744 mm1.Cross(f1, f2); //normal to the second triangle
1746 Standard_Real nn1modulus, mm1modulus;
1747 nn1modulus=sqrt(nn1.SquareModulus());
1748 mm1modulus=sqrt(mm1.SquareModulus());
1750 //-------------------------------------------------
1751 ///calculation of intersections points between triangles
1752 //-------------------------------------------------
1753 Standard_Integer NbPoints=0;
1754 IntPolyh_StartPoint SP1,SP2;
1757 if(Abs(nn1modulus)<MyConfusionPrecision) {//10.0e-20){
1761 const IntPolyh_Point n1=nn1.Divide(nn1modulus);
1763 if(NbPointsTotal<3) {
1764 IntPolyh_StartPoint Pt1,Pt2;
1765 NbPoints=TriangleEdgeContact(1,1,P1,P2,P3,e1,e2,e3,Q1,Q2,f1,n1,Pt1,Pt2);
1766 TestNbPoints(1,NbPoints,NbPointsTotal,Pt1,Pt2,SP1,SP2);
1769 if(NbPointsTotal<3) {
1770 IntPolyh_StartPoint Pt1,Pt2;
1771 NbPoints=TriangleEdgeContact(1,2,P1,P2,P3,e1,e2,e3,Q2,Q3,f2,n1,Pt1,Pt2);
1772 TestNbPoints(1,NbPoints,NbPointsTotal,Pt1,Pt2,SP1,SP2);
1775 if(NbPointsTotal<3) {
1776 IntPolyh_StartPoint Pt1,Pt2;
1777 NbPoints=TriangleEdgeContact(1,3,P1,P2,P3,e1,e2,e3,Q3,Q1,f3,n1,Pt1,Pt2);
1778 TestNbPoints(1,NbPoints,NbPointsTotal,Pt1,Pt2,SP1,SP2);
1783 if(Abs(mm1modulus)<MyConfusionPrecision) {//10.0e-20){
1787 const IntPolyh_Point m1=mm1.Divide(mm1modulus);
1789 if(NbPointsTotal<3) {
1790 IntPolyh_StartPoint Pt1,Pt2;
1791 NbPoints=TriangleEdgeContact(2,1,Q1,Q2,Q3,f1,f2,f3,P1,P2,e1,m1,Pt1,Pt2);
1792 TestNbPoints(2,NbPoints,NbPointsTotal,Pt1,Pt2,SP1,SP2);
1795 if(NbPointsTotal<3) {
1796 IntPolyh_StartPoint Pt1,Pt2;
1797 NbPoints=TriangleEdgeContact(2,2,Q1,Q2,Q3,f1,f2,f3,P2,P3,e2,m1,Pt1,Pt2);
1798 TestNbPoints(2,NbPoints,NbPointsTotal,Pt1,Pt2,SP1,SP2);
1801 if(NbPointsTotal<3) {
1802 IntPolyh_StartPoint Pt1,Pt2;
1803 NbPoints=TriangleEdgeContact(2,3,Q1,Q2,Q3,f1,f2,f3,P3,P1,e3,m1,Pt1,Pt2);
1804 TestNbPoints(2,NbPoints,NbPointsTotal,Pt1,Pt2,SP1,SP2);
1808 if (NbPointsTotal==1) {
1809 /* if( (Abs(SP1.U1()-SPInit.U1())<MyConfusionPrecision)
1810 &&(Abs(SP1.V1()-SPInit.V1())<MyConfusionPrecision) ) */
1811 if(SP1.CheckSameSP(SP2))
1819 // if ( (NbPointsTotal==2)&&( Abs(SP1.U1()-SPInit.U1())<MyConfusionPrecision)
1820 //&&( Abs(SP1.V1()-SPInit.V1())<MyConfusionPrecision) ) {
1821 if( (NbPointsTotal==2)&&(SP1.CheckSameSP(SPInit)) ) {
1822 NbPointsTotal=1;//SP1 et SPInit sont identiques
1825 // if( (NbPointsTotal==2)&&( Abs(SP2.U1()-SPInit.U1())<MyConfusionPrecision)
1826 //&&( Abs(SP2.V1()-SPInit.V1())<MyConfusionPrecision) ) {
1827 if( (NbPointsTotal==2)&&(SP2.CheckSameSP(SPInit)) ) {
1828 NbPointsTotal=1;//SP2 et SPInit sont identiques
1831 if(NbPointsTotal>1) {
1835 SPNext.SetCoupleValue(T1,T2);
1836 return (NbPointsTotal);
1838 //=======================================================================
1839 //function : NextStartingPointsResearch2
1840 //purpose : from two triangles and an intersection point I
1841 // seach the other point (if it exist).
1842 // This function is used by StartPointChain
1843 //=======================================================================
1844 Standard_Integer IntPolyh_MaillageAffinage::NextStartingPointsResearch2
1845 (const Standard_Integer T1,
1846 const Standard_Integer T2,
1847 const IntPolyh_StartPoint &SPInit,
1848 IntPolyh_StartPoint &SPNext) const
1850 Standard_Integer NbPointsTotal=0;
1851 Standard_Integer EdgeInit1=SPInit.E1();
1852 Standard_Integer EdgeInit2=SPInit.E2();
1853 if( (T1<0)||(T2<0) ) NbPointsTotal=0;
1856 const IntPolyh_Triangle &Tri1=TTriangles1[T1];
1857 const IntPolyh_Triangle &Tri2=TTriangles2[T2];
1859 const IntPolyh_Point &P1=TPoints1[Tri1.FirstPoint()];
1860 const IntPolyh_Point &P2=TPoints1[Tri1.SecondPoint()];
1861 const IntPolyh_Point &P3=TPoints1[Tri1.ThirdPoint()];
1862 const IntPolyh_Point &Q1=TPoints2[Tri2.FirstPoint()];
1863 const IntPolyh_Point &Q2=TPoints2[Tri2.SecondPoint()];
1864 const IntPolyh_Point &Q3=TPoints2[Tri2.ThirdPoint()];
1866 /* The first triangle is (p1,p2,p3). The other is (q1,q2,q3).
1867 The edges are (e1,e2,e3) and (f1,f2,f3).
1868 The normals are n1 and m1*/
1870 const IntPolyh_Point e1=P2-P1;
1871 const IntPolyh_Point e2=P3-P2;
1872 const IntPolyh_Point e3=P1-P3;
1874 const IntPolyh_Point f1=Q2-Q1;
1875 const IntPolyh_Point f2=Q3-Q2;
1876 const IntPolyh_Point f3=Q1-Q3;
1878 IntPolyh_Point nn1,mm1;
1879 nn1.Cross(e1, e2); //normal to the first triangle
1880 mm1.Cross(f1, f2); //normal to the second triangle
1882 Standard_Real nn1modulus, mm1modulus;
1883 nn1modulus=sqrt(nn1.SquareModulus());
1884 mm1modulus=sqrt(mm1.SquareModulus());
1886 //-------------------------------------------------
1887 ///calculation of intersections points between triangles
1888 //-------------------------------------------------
1890 Standard_Integer NbPoints=0;
1891 IntPolyh_StartPoint SP1,SP2;
1894 if(Abs(nn1modulus)<MyConfusionPrecision) {//10.0e-20){
1898 const IntPolyh_Point n1=nn1.Divide(nn1modulus);
1900 if( (NbPointsTotal<3)&&(EdgeInit2!=Tri2.FirstEdge()) ) {
1901 IntPolyh_StartPoint Pt1,Pt2;
1902 NbPoints=TriangleEdgeContact2(1,1,Tri1,Tri2,P1,P2,P3,e1,e2,e3,Q1,Q2,f1,n1,Pt1,Pt2);
1903 TestNbPoints(1,NbPoints,NbPointsTotal,Pt1,Pt2,SP1,SP2);
1906 if( (NbPointsTotal<3)&&(EdgeInit2!=Tri2.SecondEdge()) ) {
1907 IntPolyh_StartPoint Pt1,Pt2;
1908 NbPoints=TriangleEdgeContact2(1,2,Tri1,Tri2,P1,P2,P3,e1,e2,e3,Q2,Q3,f2,n1,Pt1,Pt2);
1909 TestNbPoints(1,NbPoints,NbPointsTotal,Pt1,Pt2,SP1,SP2);
1912 if( (NbPointsTotal<3)&&(EdgeInit2!=Tri2.ThirdEdge()) ) {
1913 IntPolyh_StartPoint Pt1,Pt2;
1914 NbPoints=TriangleEdgeContact2(1,3,Tri1,Tri2,P1,P2,P3,e1,e2,e3,Q3,Q1,f3,n1,Pt1,Pt2);
1915 TestNbPoints(1,NbPoints,NbPointsTotal,Pt1,Pt2,SP1,SP2);
1919 if(Abs(mm1modulus)<MyConfusionPrecision) {//10.0e-20){
1923 const IntPolyh_Point m1=mm1.Divide(mm1modulus);
1925 if( (NbPointsTotal<3)&&(EdgeInit1!=Tri1.FirstEdge()) ) {
1926 IntPolyh_StartPoint Pt1,Pt2;
1927 NbPoints=TriangleEdgeContact2(2,1,Tri1,Tri2,Q1,Q2,Q3,f1,f2,f3,P1,P2,e1,m1,Pt1,Pt2);
1928 TestNbPoints(2,NbPoints,NbPointsTotal,Pt1,Pt2,SP1,SP2);
1931 if( (NbPointsTotal<3)&&(EdgeInit1!=Tri1.SecondEdge()) ) {
1932 IntPolyh_StartPoint Pt1,Pt2;
1933 NbPoints=TriangleEdgeContact2(2,2,Tri1,Tri2,Q1,Q2,Q3,f1,f2,f3,P2,P3,e2,m1,Pt1,Pt2);
1934 TestNbPoints(2,NbPoints,NbPointsTotal,Pt1,Pt2,SP1,SP2);
1937 if( (NbPointsTotal<3)&&(EdgeInit1!=Tri1.ThirdEdge()) ) {
1938 IntPolyh_StartPoint Pt1,Pt2;
1939 NbPoints=TriangleEdgeContact2(2,3,Tri1,Tri2,Q1,Q2,Q3,f1,f2,f3,P3,P1,e3,m1,Pt1,Pt2);
1940 TestNbPoints(2,NbPoints,NbPointsTotal,Pt1,Pt2,SP1,SP2);
1944 if (NbPointsTotal==1) {
1945 if(SP1.CheckSameSP(SPInit))
1951 else if( (NbPointsTotal==2)&&(SP1.CheckSameSP(SPInit)) ) {
1952 NbPointsTotal=1;//SP1 et SPInit sont identiques
1955 else if( (NbPointsTotal==2)&&(SP2.CheckSameSP(SPInit)) ) {
1956 NbPointsTotal=1;//SP2 et SPInit sont identiques
1960 else if(NbPointsTotal>1) {
1964 SPNext.SetCoupleValue(T1,T2);
1965 return (NbPointsTotal);
1967 //=======================================================================
1968 //function : CalculPtsInterTriEdgeCoplanaires
1970 //=======================================================================
1971 void CalculPtsInterTriEdgeCoplanaires(const Standard_Integer TriSurfID,
1972 const IntPolyh_Point &NormaleTri,
1973 const IntPolyh_Point &PE1,
1974 const IntPolyh_Point &PE2,
1975 const IntPolyh_Point &Edge,
1976 const IntPolyh_Point &PT1,
1977 const IntPolyh_Point &PT2,
1978 const IntPolyh_Point &Cote,
1979 const Standard_Integer CoteIndex,
1980 IntPolyh_StartPoint &SP1,
1981 IntPolyh_StartPoint &SP2,
1982 Standard_Integer &NbPoints)
1984 IntPolyh_Point TestParalleles;
1985 TestParalleles.Cross(Edge,Cote);
1986 if(sqrt(TestParalleles.SquareModulus())<MyConfusionPrecision) {
1988 Per.Cross(NormaleTri,Cote);
1989 Standard_Real p1p = Per.Dot(PE1);
1990 Standard_Real p2p = Per.Dot(PE2);
1991 Standard_Real p0p = Per.Dot(PT1);
1992 if ( ( (p1p>=p0p)&&(p2p<=p0p) )||( (p1p<=p0p)&&(p2p>=p0p) ) ) {
1993 Standard_Real lambda=(p1p-p0p)/(p1p-p2p);
1994 if (lambda<-MyConfusionPrecision) {
1997 IntPolyh_Point PIE=PE1+Edge*lambda;
1998 Standard_Real alpha=RealLast();
1999 if(Cote.X()!=0) alpha=(PIE.X()-PT1.X())/Cote.X();
2000 else if (Cote.Y()!=0) alpha=(PIE.Y()-PT1.Y())/Cote.Y();
2001 else if (Cote.Z()!=0) alpha=(PIE.Z()-PT1.Z())/Cote.Z();
2005 if (alpha<-MyConfusionPrecision) {
2010 SP1.SetXYZ(PIE.X(),PIE.Y(),PIE.Z());
2012 SP1.SetUV1(PT1.U()+Cote.U()*alpha,PT1.V()+Cote.V()*alpha);
2013 SP1.SetUV2(PIE.U(),PIE.V());
2014 SP1.SetEdge1(CoteIndex);
2017 else if (TriSurfID==2) {
2018 SP1.SetUV1(PIE.U(),PIE.V());
2019 SP1.SetUV2(PT1.U()+Cote.U()*alpha,PT1.V()+Cote.V()*alpha);
2020 SP1.SetEdge2(CoteIndex);
2028 else if (NbPoints==1) {
2029 SP2.SetXYZ(PIE.X(),PIE.Y(),PIE.Z());
2031 SP2.SetUV1(PT1.U()+Cote.U()*alpha,PT1.V()+Cote.V()*alpha);
2032 SP2.SetUV2(PIE.U(),PIE.V());
2033 SP2.SetEdge1(CoteIndex);
2036 else if (TriSurfID==2) {
2037 SP2.SetUV1(PIE.U(),PIE.V());
2038 SP2.SetUV2(PT1.U()+Cote.U()*alpha,PT1.V()+Cote.V()*alpha);
2039 SP2.SetEdge2(CoteIndex);
2047 else if( (NbPoints>2)||(NbPoints<0) ) {
2053 else { //Cote et Edge paralleles, avec les rejections precedentes ils sont sur la meme droite
2054 //On projette les points sur cette droite
2055 Standard_Real pe1p= Cote.Dot(PE1);
2056 Standard_Real pe2p= Cote.Dot(PE2);
2057 Standard_Real pt1p= Cote.Dot(PT1);
2058 Standard_Real pt2p= Cote.Dot(PT2);
2060 IntPolyh_Point PEP1,PTP1,PEP2,PTP2;
2062 //PEP1 et PEP2 sont les points de contact entre le triangle et l'edge dans le repere UV de l'edge
2063 //PTP1 et PTP2 sont les correspondants respectifs a PEP1 et PEP2 dans le repere UV du triangle
2067 if ( (pt1p<pe1p)&&(pe1p<=pt2p) ) {
2069 PTP1=PT1+Cote*((pe1p-pt1p)/(pt2p-pt1p));
2073 PTP2=PT1+Cote*((pe2p-pt1p)/(pt2p-pt1p));
2077 PEP2=PE1+Edge*((pt1p-pe1p)/(pe2p-pe1p));
2082 else if( (pt2p<pe1p)&&(pe1p<=pt1p) ) {
2084 PTP1=PT1+Cote*((pt1p-pe1p)/(pt1p-pt2p));
2088 PTP2=PT1+Cote*((pe2p-pt1p)/(pt2p-pt1p));
2092 PEP2=PE1+Edge*((pt2p-pe1p)/(pe2p-pe1p));
2100 if ( (pt1p<pe2p)&&(pe2p<=pt2p) ) {
2102 PTP1=PT1+Cote*((pe2p-pt1p)/(pt2p-pt1p));
2106 PTP2=PT1+Cote*((pe1p-pt1p)/(pt2p-pt1p));
2110 PEP2=PE2+Edge*((pt1p-pe1p)/(pe2p-pe1p));
2115 else if( (pt2p<pe2p)&&(pe2p<=pt1p) ) {
2117 PTP1=PT1+Cote*((pt1p-pe2p)/(pt1p-pt2p));
2121 PTP2=PT1+Cote*((pe1p-pt1p)/(pt2p-pt1p));
2125 PEP2=PE1+Edge*((pt2p-pe1p)/(pe2p-pe1p));
2133 if (Abs(PEP1.U()-PEP2.U())<MyConfusionPrecision
2134 &&(Abs(PEP1.V()-PEP2.V())<MyConfusionPrecision) ) NbPoints=1;
2136 SP1.SetXYZ(PEP1.X(),PEP1.Y(),PEP1.Z());
2138 SP1.SetUV1(PTP1.U(),PTP1.V());
2139 SP1.SetUV2(PEP1.U(),PEP1.V());
2140 SP1.SetEdge1(CoteIndex);
2143 SP1.SetUV1(PEP1.U(),PTP1.V());
2144 SP1.SetUV2(PTP1.U(),PEP1.V());
2145 SP1.SetEdge2(CoteIndex);
2149 SP2.SetXYZ(PEP2.X(),PEP2.Y(),PEP2.Z());
2151 SP2.SetUV1(PTP2.U(),PTP2.V());
2152 SP2.SetUV2(PEP2.U(),PEP2.V());
2153 SP2.SetEdge1(CoteIndex);
2156 SP2.SetUV1(PEP2.U(),PTP2.V());
2157 SP2.SetUV2(PTP2.U(),PEP2.V());
2158 SP2.SetEdge2(CoteIndex);
2164 //=======================================================================
2165 //function : CalculPtsInterTriEdgeCoplanaires2
2167 //=======================================================================
2168 void CalculPtsInterTriEdgeCoplanaires2(const Standard_Integer TriSurfID,
2169 const IntPolyh_Point &NormaleTri,
2170 const IntPolyh_Triangle &Tri1,
2171 const IntPolyh_Triangle &Tri2,
2172 const IntPolyh_Point &PE1,
2173 const IntPolyh_Point &PE2,
2174 const IntPolyh_Point &Edge,
2175 const Standard_Integer EdgeIndex,
2176 const IntPolyh_Point &PT1,
2177 const IntPolyh_Point &PT2,
2178 const IntPolyh_Point &Cote,
2179 const Standard_Integer CoteIndex,
2180 IntPolyh_StartPoint &SP1,
2181 IntPolyh_StartPoint &SP2,
2182 Standard_Integer &NbPoints)
2184 IntPolyh_Point TestParalleles;
2185 TestParalleles.Cross(Edge,Cote);
2187 if(sqrt(TestParalleles.SquareModulus())>MyConfusionPrecision) {
2188 ///Edge and side are not parallel
2190 Per.Cross(NormaleTri,Cote);
2191 Standard_Real p1p = Per.Dot(PE1);
2192 Standard_Real p2p = Per.Dot(PE2);
2193 Standard_Real p0p = Per.Dot(PT1);
2194 ///The edge are PT1 are projected on the perpendicular of the side in the plane of the triangle
2195 if ( ( (p1p>=p0p)&&(p0p>=p2p) )||( (p1p<=p0p)&&(p0p<=p2p) ) ) {
2196 Standard_Real lambda=(p1p-p0p)/(p1p-p2p);
2197 if (lambda<-MyConfusionPrecision) {
2201 if (Abs(lambda)<MyConfusionPrecision)//lambda=0
2203 else if (Abs(lambda)>1.0-MyConfusionPrecision)//lambda=1
2206 PIE=PE1+Edge*lambda;
2208 Standard_Real alpha=RealLast();
2209 if(Cote.X()!=0) alpha=(PIE.X()-PT1.X())/Cote.X();
2210 else if (Cote.Y()!=0) alpha=(PIE.Y()-PT1.Y())/Cote.Y();
2211 else if (Cote.Z()!=0) alpha=(PIE.Z()-PT1.Z())/Cote.Z();
2216 if (alpha<-MyConfusionPrecision) {
2221 SP1.SetXYZ(PIE.X(),PIE.Y(),PIE.Z());
2223 if(Abs(alpha)<MyConfusionPrecision) {//alpha=0
2224 SP1.SetUV1(PT1.U(),PT1.V());
2225 SP1.SetUV1(PIE.U(),PIE.V());
2228 if(Abs(alpha)>1.0-MyConfusionPrecision) {//alpha=1
2229 SP1.SetUV1(PT2.U(),PT2.V());
2230 SP1.SetUV1(PIE.U(),PIE.V());
2234 SP1.SetUV1(PT1.U()+Cote.U()*alpha,PT1.V()+Cote.V()*alpha);
2235 SP1.SetUV2(PIE.U(),PIE.V());
2236 SP1.SetEdge1(Tri1.GetEdgeNumber(CoteIndex));
2237 if (Tri1.GetEdgeOrientation(CoteIndex)>0) SP1.SetLambda1(alpha);
2238 else SP1.SetLambda1(1.0-alpha);
2242 else if (TriSurfID==2) {
2243 if(Abs(alpha)<MyConfusionPrecision) {//alpha=0
2244 SP1.SetUV1(PT1.U(),PT1.V());
2245 SP1.SetUV1(PIE.U(),PIE.V());
2248 if(Abs(alpha)>1.0-MyConfusionPrecision) {//alpha=1
2249 SP1.SetUV1(PT2.U(),PT2.V());
2250 SP1.SetUV1(PIE.U(),PIE.V());
2254 SP1.SetUV1(PIE.U(),PIE.V());
2255 SP1.SetUV2(PT1.U()+Cote.U()*alpha,PT1.V()+Cote.V()*alpha);
2256 SP1.SetEdge2(Tri2.GetEdgeNumber(CoteIndex));
2257 if (Tri2.GetEdgeOrientation(CoteIndex)>0) SP1.SetLambda2(alpha);
2258 else SP1.SetLambda2(1.0-alpha);
2267 else if (NbPoints==1) {
2268 SP2.SetXYZ(PIE.X(),PIE.Y(),PIE.Z());
2270 if(Abs(alpha)<MyConfusionPrecision) {//alpha=0
2271 SP2.SetUV1(PT1.U(),PT1.V());
2272 SP2.SetUV1(PIE.U(),PIE.V());
2275 if(Abs(alpha)>1.0-MyConfusionPrecision) {//alpha=1
2276 SP2.SetUV1(PT2.U(),PT2.V());
2277 SP2.SetUV1(PIE.U(),PIE.V());
2281 SP2.SetUV1(PT1.U()+Cote.U()*alpha,PT1.V()+Cote.V()*alpha);
2282 SP2.SetUV2(PIE.U(),PIE.V());
2283 SP2.SetEdge1(Tri1.GetEdgeNumber(CoteIndex));
2284 if (Tri1.GetEdgeOrientation(CoteIndex)>0) SP2.SetLambda1(alpha);
2285 else SP2.SetLambda1(1.0-alpha);
2289 else if (TriSurfID==2) {
2290 if(Abs(alpha)<MyConfusionPrecision) {//alpha=0
2291 SP2.SetUV1(PT1.U(),PT1.V());
2292 SP2.SetUV1(PIE.U(),PIE.V());
2295 if(Abs(alpha)>1.0-MyConfusionPrecision) {//alpha=1
2296 SP2.SetUV1(PT2.U(),PT2.V());
2297 SP2.SetUV1(PIE.U(),PIE.V());
2301 SP2.SetUV1(PIE.U(),PIE.V());
2302 SP2.SetUV2(PT1.U()+Cote.U()*alpha,PT1.V()+Cote.V()*alpha);
2303 SP2.SetEdge2(Tri2.GetEdgeNumber(CoteIndex));
2304 if (Tri1.GetEdgeOrientation(CoteIndex)>0) SP2.SetLambda2(alpha);
2305 else SP2.SetLambda2(1.0-alpha);
2314 else if( (NbPoints>2)||(NbPoints<0) ) {
2321 //Side and Edge are parallel, with previous
2322 //rejections they are at the same side
2323 //The points are projected on that side
2324 Standard_Real pe1p= Cote.Dot(PE1);
2325 Standard_Real pe2p= Cote.Dot(PE2);
2326 Standard_Real pt1p= Cote.Dot(PT1);
2327 Standard_Real pt2p= Cote.Dot(PT2);
2328 Standard_Real lambda1=0., lambda2=0., alpha1=0., alpha2=0.;
2329 IntPolyh_Point PEP1,PTP1,PEP2,PTP2;
2332 if ( (pt1p<pe1p)&&(pe1p<=pt2p) ) {
2335 alpha1=((pe1p-pt1p)/(pt2p-pt1p));
2336 PTP1=PT1+Cote*alpha1;
2341 alpha2=((pe2p-pt1p)/(pt2p-pt1p));
2342 PTP2=PT1+Cote*alpha2;
2346 lambda2=((pt1p-pe1p)/(pe2p-pe1p));
2347 PEP2=PE1+Edge*lambda2;
2353 else if( (pt2p<pe1p)&&(pe1p<=pt1p) ) {
2356 alpha1=((pt1p-pe1p)/(pt1p-pt2p));
2357 PTP1=PT1+Cote*alpha1;
2362 alpha2=((pe2p-pt1p)/(pt2p-pt1p));
2363 PTP2=PT1+Cote*alpha2;
2367 lambda2=((pt2p-pe1p)/(pe2p-pe1p));
2368 PEP2=PE1+Edge*lambda2;
2377 if ( (pt1p<pe2p)&&(pe2p<=pt2p) ) {
2380 alpha1=((pe2p-pt1p)/(pt2p-pt1p));
2381 PTP1=PT1+Cote*alpha1;
2386 alpha2=((pe1p-pt1p)/(pt2p-pt1p));
2387 PTP2=PT1+Cote*alpha2;
2391 lambda2=((pt1p-pe1p)/(pe2p-pe1p));
2392 PEP2=PE2+Edge*lambda2;
2398 else if( (pt2p<pe2p)&&(pe2p<=pt1p) ) {
2401 alpha1=((pt1p-pe2p)/(pt1p-pt2p));
2402 PTP1=PT1+Cote*alpha1;
2407 alpha2=((pe1p-pt1p)/(pt2p-pt1p));
2408 PTP2=PT1+Cote*alpha2;
2412 lambda2=((pt2p-pe1p)/(pe2p-pe1p));
2413 PEP2=PE1+Edge*lambda2;
2422 SP1.SetXYZ(PEP1.X(),PEP1.Y(),PEP1.Z());
2423 if (TriSurfID==1) {///cote appartient a Tri1
2424 SP1.SetUV1(PTP1.U(),PTP1.V());
2425 SP1.SetUV2(PEP1.U(),PEP1.V());
2426 SP1.SetEdge1(Tri1.GetEdgeNumber(CoteIndex));
2428 if(Tri1.GetEdgeOrientation(CoteIndex)>0) SP1.SetLambda1(alpha1);
2429 else SP1.SetLambda1(1.0-alpha1);
2431 if(Tri2.GetEdgeOrientation(EdgeIndex)>0) SP1.SetLambda2(lambda1);
2432 else SP1.SetLambda2(1.0-lambda1);
2434 if (TriSurfID==2) {///cote appartient a Tri2
2435 SP1.SetUV1(PEP1.U(),PTP1.V());
2436 SP1.SetUV2(PTP1.U(),PEP1.V());
2437 SP1.SetEdge2(Tri2.GetEdgeNumber(CoteIndex));
2439 if(Tri2.GetEdgeOrientation(CoteIndex)>0) SP1.SetLambda1(alpha1);
2440 else SP1.SetLambda1(1.0-alpha1);
2442 if(Tri1.GetEdgeOrientation(EdgeIndex)>0) SP1.SetLambda2(lambda1);
2443 else SP1.SetLambda2(1.0-lambda1);
2446 //It is checked if PEP1!=PEP2
2447 if ( (NbPoints==2)&&(Abs(PEP1.U()-PEP2.U())<MyConfusionPrecision)
2448 &&(Abs(PEP1.V()-PEP2.V())<MyConfusionPrecision) ) NbPoints=1;
2450 SP2.SetXYZ(PEP2.X(),PEP2.Y(),PEP2.Z());
2452 SP2.SetUV1(PTP2.U(),PTP2.V());
2453 SP2.SetUV2(PEP2.U(),PEP2.V());
2454 SP2.SetEdge1(Tri1.GetEdgeNumber(CoteIndex));
2456 if(Tri1.GetEdgeOrientation(CoteIndex)>0) SP2.SetLambda1(alpha1);
2457 else SP2.SetLambda1(1.0-alpha1);
2459 if(Tri2.GetEdgeOrientation(EdgeIndex)>0) SP2.SetLambda2(lambda1);
2460 else SP2.SetLambda2(1.0-lambda1);
2463 SP2.SetUV1(PEP2.U(),PTP2.V());
2464 SP2.SetUV2(PTP2.U(),PEP2.V());
2465 SP2.SetEdge2(Tri2.GetEdgeNumber(CoteIndex));
2467 if(Tri1.GetEdgeOrientation(CoteIndex)>0) SP2.SetLambda1(alpha1);
2468 else SP2.SetLambda1(1.0-alpha1);
2470 if(Tri2.GetEdgeOrientation(EdgeIndex)>0) SP2.SetLambda2(lambda1);
2471 else SP2.SetLambda2(1.0-lambda1);
2476 //Filter if the point is placed on top, the edge is set to -1
2478 if(Abs(SP1.Lambda1())<MyConfusionPrecision)
2480 if(Abs(SP1.Lambda1()-1)<MyConfusionPrecision)
2482 if(Abs(SP1.Lambda2())<MyConfusionPrecision)
2484 if(Abs(SP1.Lambda2()-1)<MyConfusionPrecision)
2488 if(Abs(SP2.Lambda1())<MyConfusionPrecision)
2490 if(Abs(SP2.Lambda1()-1)<MyConfusionPrecision)
2492 if(Abs(SP2.Lambda2())<MyConfusionPrecision)
2494 if(Abs(SP2.Lambda2()-1)<MyConfusionPrecision)
2498 //=======================================================================
2499 //function : TriangleEdgeContact
2501 //=======================================================================
2502 Standard_Integer IntPolyh_MaillageAffinage::TriangleEdgeContact
2503 (const Standard_Integer TriSurfID,
2504 const Standard_Integer EdgeIndex,
2505 const IntPolyh_Point &PT1,
2506 const IntPolyh_Point &PT2,
2507 const IntPolyh_Point &PT3,
2508 const IntPolyh_Point &Cote12,
2509 const IntPolyh_Point &Cote23,
2510 const IntPolyh_Point &Cote31,
2511 const IntPolyh_Point &PE1,const IntPolyh_Point &PE2,
2512 const IntPolyh_Point &Edge,
2513 const IntPolyh_Point &NormaleT,
2514 IntPolyh_StartPoint &SP1,IntPolyh_StartPoint &SP2) const
2517 Standard_Real lambda =0.;
2518 Standard_Real alpha =0.;
2519 Standard_Real beta =0.;
2522 //The edge, on which the point is located, is known.
2524 SP1.SetEdge2(EdgeIndex);
2525 SP2.SetEdge2(EdgeIndex);
2527 else if (TriSurfID==2) {
2528 SP1.SetEdge1(EdgeIndex);
2529 SP2.SetEdge1(EdgeIndex);
2535 /**The edge is projected on the normal of the triangle if yes
2536 --> free intersection (point I)--> start point is found*/
2537 Standard_Integer NbPoints=0;
2539 if(NormaleT.SquareModulus()==0) {
2542 else if( (Cote12.SquareModulus()==0)
2543 ||(Cote23.SquareModulus()==0)
2544 ||(Cote31.SquareModulus()==0) ) {
2547 else if(Edge.SquareModulus()==0) {
2551 Standard_Real pe1 = NormaleT.Dot(PE1);
2552 Standard_Real pe2 = NormaleT.Dot(PE2);
2553 Standard_Real pt1 = NormaleT.Dot(PT1);
2555 // PE1I = lambda.Edge
2557 if( (Abs(pe1-pe2)<MyConfusionPrecision)&&(Abs(pe1-pt1)<MyConfusionPrecision) ) {
2558 //edge and triangle are coplanar (two contact points maximum)
2560 //The tops of the triangle are projected on the perpendicular of the edge
2562 IntPolyh_Point PerpEdge;
2563 PerpEdge.Cross(NormaleT,Edge);
2564 Standard_Real pp1 = PerpEdge.Dot(PT1);
2565 Standard_Real pp2 = PerpEdge.Dot(PT2);
2566 Standard_Real pp3 = PerpEdge.Dot(PT3);
2567 Standard_Real ppe1 = PerpEdge.Dot(PE1);
2569 if ( ( (pp1>ppe1)&&(pp2<=ppe1)&&(pp3<=ppe1) ) || ( (pp1<ppe1)&&(pp2>=ppe1)&&(pp3>=ppe1) ) ){
2570 //there are two sides (commun top PT1) that can cut the edge
2573 CalculPtsInterTriEdgeCoplanaires(TriSurfID,NormaleT,
2574 PE1,PE2,Edge,PT1,PT2,
2575 Cote12,1,SP1,SP2,NbPoints);
2577 if ( (NbPoints>1)&&(Abs(SP2.U1()-SP1.U1())<MyConfusionPrecision)
2578 &&(Abs(SP1.V1()-SP2.V1())<MyConfusionPrecision) ) NbPoints=1;
2582 CalculPtsInterTriEdgeCoplanaires(TriSurfID,NormaleT,
2583 PE1,PE2,Edge,PT3,PT1,
2584 Cote31,3,SP1,SP2,NbPoints);
2588 if ( (NbPoints>1)&&(Abs(SP1.U1()-SP2.U1())<MyConfusionPrecision)
2589 &&(Abs(SP1.V2()-SP2.V1())<MyConfusionPrecision) ) NbPoints=1;
2590 if (NbPoints>=2) return(NbPoints);
2592 else if ( ( ( (pp2>ppe1)&&(pp1<=ppe1)&&(pp3<=ppe1) ) || ( (pp2<ppe1)&&(pp1>=ppe1)&&(pp3>=ppe1) ) )
2594 //there are two sides (common top PT2) that can cut the edge
2597 CalculPtsInterTriEdgeCoplanaires(TriSurfID,NormaleT,
2598 PE1,PE2,Edge,PT1,PT2,
2599 Cote12,1,SP1,SP2,NbPoints);
2601 if ( (NbPoints>1)&&(Abs(SP2.U1()-SP1.U1())<MyConfusionPrecision)
2602 &&(Abs(SP1.V1()-SP2.V1())<MyConfusionPrecision) ) NbPoints=1;
2605 if(NbPoints<2) CalculPtsInterTriEdgeCoplanaires(TriSurfID,NormaleT,
2606 PE1,PE2,Edge,PT2,PT3,
2607 Cote23,2,SP1,SP2,NbPoints);
2609 if ( (NbPoints>1)&&(Abs(SP2.U1()-SP1.U1())<MyConfusionPrecision)
2610 &&(Abs(SP1.V1()-SP2.V1())<MyConfusionPrecision) ) NbPoints=1;
2611 if (NbPoints>=2) return(NbPoints);
2613 else if ( ( (pp3>ppe1)&&(pp1<=ppe1)&&(pp2<=ppe1) ) || ( (pp3<ppe1)&&(pp1>=ppe1)&&(pp2>=ppe1) )
2615 //there are two sides (common top PT3) that can cut the edge
2618 CalculPtsInterTriEdgeCoplanaires(TriSurfID,NormaleT,
2619 PE1,PE2,Edge,PT3,PT1,Cote31,
2620 3,SP1,SP2,NbPoints);
2622 if ( (NbPoints>1)&&(Abs(SP2.U1()-SP1.U1())<MyConfusionPrecision)
2623 &&(Abs(SP1.V1()-SP2.V1())<MyConfusionPrecision) ) NbPoints=1;
2626 if (NbPoints<2) CalculPtsInterTriEdgeCoplanaires(TriSurfID,NormaleT,
2627 PE1,PE2,Edge,PT2,PT3,
2628 Cote23,2,SP1,SP2,NbPoints);
2630 if ( (NbPoints>1)&&(Abs(SP2.U1()-SP1.U1())<MyConfusionPrecision)
2631 &&(Abs(SP2.V1()-SP1.V1())<MyConfusionPrecision) ) NbPoints=1;
2632 if (NbPoints>=2) return(NbPoints);
2636 //------------------------------------------------------
2637 // edge and triangle NON COPLANAR (a contact point)
2638 //------------------------------------------------------
2639 else if( ( (pe1>=pt1)&&(pe2<=pt1) ) || ( (pe1<=pt1)&&(pe2>=pt1) ) ) {
2640 lambda=(pe1-pt1)/(pe1-pe2);
2642 if (lambda<-MyConfusionPrecision) {
2645 else if (Abs(lambda)<MyConfusionPrecision) {//lambda==0
2647 if(TriSurfID==1) SP1.SetEdge2(0);
2648 else SP1.SetEdge1(0);
2650 else if (Abs(lambda-1.0)<MyConfusionPrecision) {//lambda==1
2652 if(TriSurfID==1) SP1.SetEdge2(0);
2653 else SP1.SetEdge1(0);
2657 if(TriSurfID==1) SP1.SetEdge2(EdgeIndex);
2658 else SP1.SetEdge1(EdgeIndex);
2661 if(Abs(Cote23.X())>MyConfusionPrecision) {
2662 Standard_Real D=(Cote12.Y()-Cote12.X()*Cote23.Y()/Cote23.X());
2663 if(D!=0) alpha = (PI.Y()-PT1.Y()-(PI.X()-PT1.X())*Cote23.Y()/Cote23.X())/D;
2667 if ((alpha<-MyConfusionPrecision)||(alpha>(1.0+MyConfusionPrecision))) return(0);
2668 else beta = (PI.X()-PT1.X()-alpha*Cote12.X())/Cote23.X();
2671 else if (Abs(Cote12.X())>MyConfusionPrecision) { //On a Cote23.X()==0
2672 alpha = (PI.X()-PT1.X())/Cote12.X();
2674 if ((alpha<-MyConfusionPrecision)||(alpha>(1.0+MyConfusionPrecision))) return(0);
2676 else if (Abs(Cote23.Y())>MyConfusionPrecision) beta = (PI.Y()-PT1.Y()-alpha*Cote12.Y())/Cote23.Y();
2677 else if (Abs(Cote23.Z())>MyConfusionPrecision) beta = (PI.Z()-PT1.Z()-alpha*Cote12.Z())/Cote23.Z();
2683 else if (Abs(Cote23.Y())>MyConfusionPrecision) {
2684 //On a Cote23.X()==0 et Cote12.X()==0 ==> equation can't be used
2685 Standard_Real D=(Cote12.Z()-Cote12.Y()*Cote23.Z()/Cote23.Y());
2687 if(D!=0) alpha = (PI.Z()-PT1.Z()-(PI.Y()-PT1.Y())*Cote23.Z()/Cote23.Y())/D;
2692 if ( (alpha<-MyConfusionPrecision)||(alpha>(1.0+MyConfusionPrecision)) ) return(0);
2693 else beta = (PI.Y()-PT1.Y()-alpha*Cote12.Y())/Cote23.Y();
2696 else if (Abs(Cote12.Y())>MyConfusionPrecision) {
2697 //On a Cote23.X()==0, Cote12.X()==0 et Cote23.Y()==0
2698 alpha = (PI.Y()-PT1.Y())/Cote12.Y();
2700 if ((Abs(alpha)<MyConfusionPrecision)||(Abs(alpha-1.0)<MyConfusionPrecision)) return(0);
2702 else if (Abs(Cote23.Z())>MyConfusionPrecision) beta = (PI.Z()-PT1.Z()-alpha*Cote12.Z())/Cote23.Z();
2709 else { //two equations of three can't be used
2715 if( (beta<-MyConfusionPrecision)||(beta>(alpha+MyConfusionPrecision)) ) return(0);
2717 SP1.SetXYZ(PI.X(),PI.Y(),PI.Z());
2720 SP1.SetUV2(PI.U(),PI.V());
2721 SP1.SetUV1(PT1.U()+Cote12.U()*alpha+Cote23.U()*beta, PT1.V()+Cote12.V()*alpha+Cote23.V()*beta);
2723 if (beta<MyConfusionPrecision) {//beta==0 && alpha
2725 SP1.SetLambda1(alpha);
2727 if (Abs(beta-alpha)<MyConfusionPrecision) {//beta==alpha
2729 SP1.SetLambda1(1.0-alpha);
2731 if (Abs(alpha-1)<MyConfusionPrecision)//alpha==1
2733 if (alpha<MyConfusionPrecision) {//alpha=0 --> beta==0
2734 SP1.SetXYZ(PT1.X(),PT1.Y(),PT1.Z());
2735 SP1.SetUV1(PT1.U(),PT1.V());
2738 if ( (beta<MyConfusionPrecision)&&(Abs(1-alpha)<MyConfusionPrecision) ) {//beta==0 alpha==1
2739 SP1.SetXYZ(PT2.X(),PT2.Y(),PT2.Z());
2740 SP1.SetUV1(PT2.U(),PT2.V());
2743 if ( (Abs(beta-1)<MyConfusionPrecision)||(Abs(1-alpha)<MyConfusionPrecision) ) {//beta==1 alpha==1
2744 SP1.SetXYZ(PT3.X(),PT3.Y(),PT3.Z());
2745 SP1.SetUV1(PT3.U(),PT3.V());
2749 else if(TriSurfID==2) {
2750 SP1.SetUV1(PI.U(),PI.V());
2751 SP1.SetUV2(PT1.U()+Cote12.U()*alpha+Cote23.U()*beta, PT1.V()+Cote12.V()*alpha+Cote23.V()*beta);
2753 if (beta<MyConfusionPrecision) { //beta==0
2756 if (Abs(beta-alpha)<MyConfusionPrecision)//beta==alpha
2758 if (Abs(alpha-1)<MyConfusionPrecision)//alpha==1
2760 if (alpha<MyConfusionPrecision) {//alpha=0 --> beta==0
2761 SP1.SetXYZ(PT1.X(),PT1.Y(),PT1.Z());
2762 SP1.SetUV2(PT1.U(),PT1.V());
2765 if ( (beta<MyConfusionPrecision)&&(Abs(1-alpha)<MyConfusionPrecision) ) {//beta==0 alpha==1
2766 SP1.SetXYZ(PT2.X(),PT2.Y(),PT2.Z());
2767 SP1.SetUV2(PT2.U(),PT2.V());
2770 if ( (Abs(beta-1)<MyConfusionPrecision)||(Abs(1-alpha)<MyConfusionPrecision) ) {//beta==1 alpha==1
2771 SP1.SetXYZ(PT3.X(),PT3.Y(),PT3.Z());
2772 SP1.SetUV2(PT3.U(),PT3.V());
2786 //=======================================================================
2787 //function : TriangleEdgeContact2
2789 //=======================================================================
2790 Standard_Integer IntPolyh_MaillageAffinage::TriangleEdgeContact2
2791 (const Standard_Integer TriSurfID,
2792 const Standard_Integer EdgeIndex,
2793 const IntPolyh_Triangle &Tri1,
2794 const IntPolyh_Triangle &Tri2,
2795 const IntPolyh_Point &PT1,
2796 const IntPolyh_Point &PT2,
2797 const IntPolyh_Point &PT3,
2798 const IntPolyh_Point &Cote12,
2799 const IntPolyh_Point &Cote23,
2800 const IntPolyh_Point &Cote31,
2801 const IntPolyh_Point &PE1,const IntPolyh_Point &PE2,
2802 const IntPolyh_Point &Edge,
2803 const IntPolyh_Point &NormaleT,
2804 IntPolyh_StartPoint &SP1,IntPolyh_StartPoint &SP2) const
2807 Standard_Real lambda =0., alpha =0., beta =0.;
2809 //One of edges on which the point is located is known
2811 SP1.SetEdge2(Tri2.GetEdgeNumber(EdgeIndex));
2812 SP2.SetEdge2(Tri2.GetEdgeNumber(EdgeIndex));
2814 else if (TriSurfID==2) {
2815 SP1.SetEdge1(Tri1.GetEdgeNumber(EdgeIndex));
2816 SP2.SetEdge1(Tri1.GetEdgeNumber(EdgeIndex));
2822 /**The edge is projected on the normal in the triangle if yes
2823 --> a free intersection (point I)--> Start point is found */
2824 Standard_Integer NbPoints=0;
2825 if(NormaleT.SquareModulus()==0) {
2828 else if( (Cote12.SquareModulus()==0)
2829 ||(Cote23.SquareModulus()==0)
2830 ||(Cote31.SquareModulus()==0) ) {
2833 else if(Edge.SquareModulus()==0) {
2837 Standard_Real pe1 = NormaleT.Dot(PE1);
2838 Standard_Real pe2 = NormaleT.Dot(PE2);
2839 Standard_Real pt1 = NormaleT.Dot(PT1);
2841 // PE1I = lambda.Edge
2843 if( (Abs(pe1-pt1)<MyConfusionPrecision)&&(Abs(pe2-pt1)<MyConfusionPrecision)) {
2844 //edge and triangle are coplanar (two contact points at maximum)
2847 //the tops of the triangle are projected on the perpendicular to the edge
2848 IntPolyh_Point PerpEdge;
2849 PerpEdge.Cross(NormaleT,Edge);
2850 Standard_Real pp1 = PerpEdge.Dot(PT1);
2851 Standard_Real pp2 = PerpEdge.Dot(PT2);
2852 Standard_Real pp3 = PerpEdge.Dot(PT3);
2853 Standard_Real ppe1 = PerpEdge.Dot(PE1);
2856 if ( (Abs(pp1-pp2)<MyConfusionPrecision)&&(Abs(pp1-pp3)<MyConfusionPrecision) ) {
2860 if ( ( (pp1>=ppe1)&&(pp2<=ppe1)&&(pp3<=ppe1) ) || ( (pp1<=ppe1)&&(pp2>=ppe1)&&(pp3>=ppe1) ) ){
2861 //there are two sides (common top PT1) that can cut the edge
2864 CalculPtsInterTriEdgeCoplanaires2(TriSurfID,NormaleT,Tri1,Tri2,PE1,PE2,Edge,EdgeIndex,
2865 PT1,PT2,Cote12,1,SP1,SP2,NbPoints);
2867 if ( (NbPoints>1)&&(Abs(SP2.U1()-SP1.U1())<MyConfusionPrecision)
2868 &&(Abs(SP1.V1()-SP2.V1())<MyConfusionPrecision) ) NbPoints=1;
2871 if (NbPoints<2) CalculPtsInterTriEdgeCoplanaires2(TriSurfID,NormaleT,Tri1,Tri2,PE1,PE2,Edge,EdgeIndex,
2872 PT3,PT1,Cote31,3,SP1,SP2,NbPoints);
2875 if ( (NbPoints>1)&&(Abs(SP1.U1()-SP2.U1())<MyConfusionPrecision)
2876 &&(Abs(SP1.V2()-SP2.V1())<MyConfusionPrecision) ) NbPoints=1;
2877 if (NbPoints>=2) return(NbPoints);
2879 else if ( ( ( (pp2>=ppe1)&&(pp1<=ppe1)&&(pp3<=ppe1) ) || ( (pp2<=ppe1)&&(pp1>=ppe1)&&(pp3>=ppe1) ) )
2881 //there are two sides (common top PT2) that can cut the edge
2884 CalculPtsInterTriEdgeCoplanaires2(TriSurfID,NormaleT,Tri1,Tri2,PE1,PE2,Edge,EdgeIndex,
2885 PT1,PT2,Cote12,1,SP1,SP2,NbPoints);
2887 if ( (NbPoints>1)&&(Abs(SP2.U1()-SP1.U1())<MyConfusionPrecision)
2888 &&(Abs(SP1.V1()-SP2.V1())<MyConfusionPrecision) ) NbPoints=1;
2891 if(NbPoints<2) CalculPtsInterTriEdgeCoplanaires2(TriSurfID,NormaleT,Tri1,Tri2,PE1,PE2,Edge,EdgeIndex,
2892 PT2,PT3,Cote23,2,SP1,SP2,NbPoints);
2894 if ( (NbPoints>1)&&(Abs(SP2.U1()-SP1.U1())<MyConfusionPrecision)
2895 &&(Abs(SP1.V1()-SP2.V1())<MyConfusionPrecision) ) NbPoints=1;
2896 if (NbPoints>=2) return(NbPoints);
2898 else if ( ( (pp3>=ppe1)&&(pp1<=ppe1)&&(pp2<=ppe1) ) || ( (pp3<=ppe1)&&(pp1>=ppe1)&&(pp2>=ppe1) )
2900 //there are two sides (common top PT3) that can cut the edge
2903 CalculPtsInterTriEdgeCoplanaires2(TriSurfID,NormaleT,Tri1,Tri2,PE1,PE2,Edge,EdgeIndex,
2904 PT3,PT1,Cote31,3,SP1,SP2,NbPoints);
2906 if ( (NbPoints>1)&&(Abs(SP2.U1()-SP1.U1())<MyConfusionPrecision)
2907 &&(Abs(SP1.V1()-SP2.V1())<MyConfusionPrecision) ) NbPoints=1;
2910 if (NbPoints<2) CalculPtsInterTriEdgeCoplanaires2(TriSurfID,NormaleT,Tri1,Tri2,PE1,PE2,Edge,EdgeIndex,
2911 PT2,PT3,Cote23,2,SP1,SP2,NbPoints);
2913 if ( (NbPoints>1)&&(Abs(SP2.U1()-SP1.U1())<MyConfusionPrecision)
2914 &&(Abs(SP2.V1()-SP1.V1())<MyConfusionPrecision) ) NbPoints=1;
2915 if (NbPoints>=2) return(NbPoints);
2919 //------------------------------------------------------
2920 // NON COPLANAR edge and triangle (a contact point)
2921 //------------------------------------------------------
2922 else if( ( (pe1>=pt1)&&(pt1>=pe2) ) || ( (pe1<=pt1)&&(pt1<=pe2) ) ) { //
2923 lambda=(pe1-pt1)/(pe1-pe2);
2925 if (lambda<-MyConfusionPrecision) {
2928 else if (Abs(lambda)<MyConfusionPrecision) {//lambda==0
2930 if(TriSurfID==1) SP1.SetEdge2(-1);
2931 else SP1.SetEdge1(-1);
2933 else if (Abs(lambda-1.0)<MyConfusionPrecision) {//lambda==1
2935 if(TriSurfID==1) SP1.SetEdge2(-1);
2936 else SP1.SetEdge1(-1);
2941 if(Tri2.GetEdgeOrientation(EdgeIndex)>0)
2942 SP1.SetLambda2(lambda);
2943 else SP1.SetLambda2(1.0-lambda);
2945 if(Tri1.GetEdgeOrientation(EdgeIndex)>0)
2946 SP1.SetLambda1(lambda);
2947 else SP1.SetLambda1(1.0-lambda);
2951 Standard_Real Cote23X=Cote23.X();
2952 Standard_Real D1=0.0;
2953 Standard_Real D3,D4;
2955 //Combination Eq1 Eq2
2956 if(Abs(Cote23X)>MyConfusionPrecision) {
2957 D1=Cote12.Y()-Cote12.X()*Cote23.Y()/Cote23X;
2959 if(Abs(D1)>MyConfusionPrecision) {
2960 alpha = ( PI.Y()-PT1.Y()-(PI.X()-PT1.X())*Cote23.Y()/Cote23X )/D1;
2962 ///It is checked if 1.0>=alpha>=0.0
2963 if ((alpha<-MyConfusionPrecision)||(alpha>(1.0+MyConfusionPrecision))) return(0);
2964 else beta = (PI.X()-PT1.X()-alpha*Cote12.X())/Cote23X;
2966 //Combination Eq1 and Eq2 with Cote23.X()==0
2967 else if ( (Abs(Cote12.X())>MyConfusionPrecision)
2968 &&(Abs(Cote23X)<MyConfusionPrecision) ) { //There is Cote23.X()==0
2969 alpha = (PI.X()-PT1.X())/Cote12.X();
2971 if ((alpha<-MyConfusionPrecision)||(alpha>(1.0+MyConfusionPrecision))) return(0);
2973 else if (Abs(Cote23.Y())>MyConfusionPrecision) beta = (PI.Y()-PT1.Y()-alpha*Cote12.Y())/Cote23.Y();
2974 else if (Abs(Cote23.Z())>MyConfusionPrecision) beta = (PI.Z()-PT1.Z()-alpha*Cote12.Z())/Cote23.Z();
2979 //Combination Eq1 and Eq3
2980 else if ( (Abs(Cote23.X())>MyConfusionPrecision)
2981 &&(Abs( D3= (Cote12.Z()-Cote12.X()*Cote23.Z()/Cote23.X()) ) > MyConfusionPrecision) ) {
2983 alpha = (PI.Z()-PT1.Z()-(PI.X()-PT1.X())*Cote23.Z()/Cote23.X())/D3;
2985 if ( (alpha<-MyConfusionPrecision)||(alpha>(1.0+MyConfusionPrecision)) ) return(0);
2986 else beta = (PI.X()-PT1.X()-alpha*Cote12.X())/Cote23.X();
2988 //Combination Eq2 and Eq3
2989 else if ( (Abs(Cote23.Y())>MyConfusionPrecision)
2990 &&(Abs( D4= (Cote12.Z()-Cote12.Y()*Cote23.Z()/Cote23.Y()) ) > MyConfusionPrecision) ) {
2992 alpha = (PI.Z()-PT1.Z()-(PI.Y()-PT1.Y())*Cote23.Z()/Cote23.Y())/D4;
2994 if ( (alpha<-MyConfusionPrecision)||(alpha>(1.0+MyConfusionPrecision)) ) return(0);
2995 else beta = (PI.Y()-PT1.Y()-alpha*Cote12.Y())/Cote23.Y();
2997 //Combination Eq2 and Eq3 with Cote23.Y()==0
2998 else if ( (Abs(Cote12.Y())>MyConfusionPrecision)
2999 && (Abs(Cote23.Y())<MyConfusionPrecision) ) {
3000 alpha = (PI.Y()-PT1.Y())/Cote12.Y();
3002 if ( (alpha<-MyConfusionPrecision)||(alpha>(1.0+MyConfusionPrecision)) ) return(0);
3004 else if (Abs(Cote23.Z())>MyConfusionPrecision) beta = (PI.Z()-PT1.Z()-alpha*Cote12.Z())/Cote23.Z();
3007 printf("\nCote PT2PT3 nul1\n");
3012 //Combination Eq1 and Eq3 with Cote23.Z()==0
3013 else if ( (Abs(Cote12.Z())>MyConfusionPrecision)
3014 && (Abs(Cote23.Z())<MyConfusionPrecision) ) {
3015 alpha = (PI.Z()-PT1.Z())/Cote12.Z();
3017 if ( (alpha<-MyConfusionPrecision)||(alpha>(1.0+MyConfusionPrecision)) ) return(0);
3019 else if (Abs(Cote23.X())>MyConfusionPrecision) beta = (PI.X()-PT1.X()-alpha*Cote12.X())/Cote23.X();
3026 else { //Particular case not processed ?
3032 if( (beta<-MyConfusionPrecision)||(beta>(alpha+MyConfusionPrecision)) ) return(0);
3034 SP1.SetXYZ(PI.X(),PI.Y(),PI.Z());
3037 SP1.SetUV2(PI.U(),PI.V());
3038 SP1.SetUV1(PT1.U()+Cote12.U()*alpha+Cote23.U()*beta, PT1.V()+Cote12.V()*alpha+Cote23.V()*beta);
3040 if (alpha<MyConfusionPrecision) {//alpha=0 --> beta==0
3041 SP1.SetXYZ(PT1.X(),PT1.Y(),PT1.Z());
3042 SP1.SetUV1(PT1.U(),PT1.V());
3045 else if ( (beta<MyConfusionPrecision)&&(Abs(1-alpha)<MyConfusionPrecision) ) {//beta==0 alpha==1
3046 SP1.SetXYZ(PT2.X(),PT2.Y(),PT2.Z());
3047 SP1.SetUV1(PT2.U(),PT2.V());
3050 else if ( (Abs(beta-1)<MyConfusionPrecision)&&(Abs(1-alpha)<MyConfusionPrecision) ) {//beta==1 alpha==1
3051 SP1.SetXYZ(PT3.X(),PT3.Y(),PT3.Z());
3052 SP1.SetUV1(PT3.U(),PT3.V());
3055 else if (beta<MyConfusionPrecision) {//beta==0
3056 SP1.SetEdge1(Tri1.GetEdgeNumber(1));
3057 if(Tri1.GetEdgeOrientation(1)>0)
3058 SP1.SetLambda1(alpha);
3059 else SP1.SetLambda1(1.0-alpha);
3061 else if (Abs(beta-alpha)<MyConfusionPrecision) {//beta==alpha
3062 SP1.SetEdge1(Tri1.GetEdgeNumber(3));
3063 if(Tri1.GetEdgeOrientation(3)>0)
3064 SP1.SetLambda1(1.0-alpha);
3065 else SP1.SetLambda1(alpha);
3067 else if (Abs(alpha-1)<MyConfusionPrecision) {//alpha==1
3068 SP1.SetEdge1(Tri1.GetEdgeNumber(2));
3069 if(Tri1.GetEdgeOrientation(2)>0)
3070 SP1.SetLambda1(beta);
3071 else SP1.SetLambda1(1.0-beta);
3074 else if(TriSurfID==2) {
3075 SP1.SetUV1(PI.U(),PI.V());
3076 SP1.SetUV2(PT1.U()+Cote12.U()*alpha+Cote23.U()*beta, PT1.V()+Cote12.V()*alpha+Cote23.V()*beta);
3078 if (alpha<MyConfusionPrecision) {//alpha=0 --> beta==0
3079 SP1.SetXYZ(PT1.X(),PT1.Y(),PT1.Z());
3080 SP1.SetUV2(PT1.U(),PT1.V());
3083 else if ( (beta<MyConfusionPrecision)&&(Abs(1-alpha)<MyConfusionPrecision) ) {//beta==0 alpha==1
3084 SP1.SetXYZ(PT2.X(),PT2.Y(),PT2.Z());
3085 SP1.SetUV2(PT2.U(),PT2.V());
3088 else if ( (Abs(beta-1)<MyConfusionPrecision)&&(Abs(1-alpha)<MyConfusionPrecision) ) {//beta==1 alpha==1
3089 SP1.SetXYZ(PT3.X(),PT3.Y(),PT3.Z());
3090 SP1.SetUV2(PT3.U(),PT3.V());
3093 else if (beta<MyConfusionPrecision) { //beta==0
3094 SP1.SetEdge2(Tri2.GetEdgeNumber(1));
3095 if(Tri2.GetEdgeOrientation(1)>0)
3096 SP1.SetLambda2(alpha);
3097 else SP1.SetLambda2(1.0-alpha);
3099 else if (Abs(beta-alpha)<MyConfusionPrecision) {//beta==alpha
3100 SP1.SetEdge2(Tri2.GetEdgeNumber(3));
3101 if(Tri2.GetEdgeOrientation(3)>0)
3102 SP1.SetLambda2(1.0-alpha);
3103 else SP1.SetLambda2(alpha);
3105 else if (Abs(alpha-1)<MyConfusionPrecision) {//alpha==1
3106 SP1.SetEdge2(Tri2.GetEdgeNumber(2));
3107 if(Tri2.GetEdgeOrientation(2)>0)
3108 SP1.SetLambda2(alpha);
3109 else SP1.SetLambda2(1.0-alpha);
3121 //=======================================================================
3122 //function : TriangleComparePSP
3123 //purpose : The same as TriangleCompare, plus compute the
3124 // StartPoints without chaining them.
3125 //=======================================================================
3126 Standard_Integer IntPolyh_MaillageAffinage::TriangleComparePSP ()
3128 Standard_Integer CpteurTab=0;
3129 Standard_Integer CpteurTabSP=0;
3130 Standard_Real CoupleAngle=-2.0;
3131 const Standard_Integer FinTT1 = TTriangles1.NbItems();
3132 const Standard_Integer FinTT2 = TTriangles2.NbItems();
3134 for(Standard_Integer i_S1=0; i_S1<FinTT1; i_S1++) {
3135 IntPolyh_Triangle &Triangle1 = TTriangles1[i_S1];
3136 if ((Triangle1.IndiceIntersectionPossible() == 0) ||
3137 (Triangle1.GetFleche() < 0.))
3139 for(Standard_Integer i_S2=0; i_S2<FinTT2; i_S2++){
3140 IntPolyh_Triangle &Triangle2 = TTriangles2[i_S2];
3141 if ((Triangle2.IndiceIntersectionPossible() != 0) &&
3142 (Triangle2.GetFleche() >= 0.)) {
3143 IntPolyh_StartPoint SP1, SP2;
3144 //If a triangle is dead or not in BSB, comparison is not possible
3146 Standard_Integer iDeg1, iDeg2, iDeg3, iDeg;
3148 const IntPolyh_Point& P1=TPoints1[Triangle1.FirstPoint()];
3149 const IntPolyh_Point& P2=TPoints1[Triangle1.SecondPoint()];
3150 const IntPolyh_Point& P3=TPoints1[Triangle1.ThirdPoint()];
3151 iDeg1=(P1.Degenerated()) ? 1 : 0;
3152 iDeg2=(P2.Degenerated()) ? 1 : 0;
3153 iDeg3=(P3.Degenerated()) ? 1 : 0;
3154 iDeg=iDeg1+iDeg2+iDeg3;
3159 const IntPolyh_Point& Q1=TPoints2[Triangle2.FirstPoint()];
3160 const IntPolyh_Point& Q2=TPoints2[Triangle2.SecondPoint()];
3161 const IntPolyh_Point& Q3=TPoints2[Triangle2.ThirdPoint()];
3162 iDeg1=(Q1.Degenerated()) ? 1 : 0;
3163 iDeg2=(Q2.Degenerated()) ? 1 : 0;
3164 iDeg3=(Q3.Degenerated()) ? 1 : 0;
3165 iDeg=iDeg1+iDeg2+iDeg3;
3170 if (TriContact(P1, P2, P3, Q1, Q2, Q3, CoupleAngle)) {
3171 Triangle1.SetIndiceIntersection(1);//The triangle is cut by another
3172 Triangle2.SetIndiceIntersection(1);
3174 Standard_Integer NbPoints;
3175 NbPoints=StartingPointsResearch(i_S1,i_S2,SP1, SP2);
3181 if ( (NbPoints>0)&&(NbPoints<3) ) {
3182 SP1.SetCoupleValue(i_S1,i_S2);
3183 TStartPoints[CpteurTabSP]=SP1;
3190 SP2.SetCoupleValue(i_S1,i_S2);
3191 TStartPoints[CpteurTabSP]=SP2;
3205 return(CpteurTabSP);
3207 //=======================================================================
3208 //function : TriangleCompare
3209 //purpose : Analyze each couple of triangles from the two --
3210 // array of triangles, to see if they are in
3211 // contact, and compute the incidence. Then put
3212 // couples in contact in the array of couples
3213 //=======================================================================
3214 Standard_Integer IntPolyh_MaillageAffinage::TriangleCompare ()
3216 Standard_Integer CpteurTab=0;
3218 const Standard_Integer FinTT1 = TTriangles1.NbItems();
3219 const Standard_Integer FinTT2 = TTriangles2.NbItems();
3221 Standard_Integer TTClimit = 200;
3222 Standard_Integer NbTTC = FinTT1 * FinTT2 / 10;
3223 if (NbTTC < TTClimit)
3225 TTrianglesContacts.Init(NbTTC);
3227 //TTrianglesContacts.Init(FinTT1 * FinTT2 / 10);
3229 Standard_Real CoupleAngle=-2.0;
3230 for(Standard_Integer i_S1=0; i_S1<FinTT1; i_S1++) {
3231 IntPolyh_Triangle &Triangle1 = TTriangles1[i_S1];
3232 if ((Triangle1.IndiceIntersectionPossible() == 0) ||
3233 (Triangle1.GetFleche() < 0.))
3235 for(Standard_Integer i_S2=0; i_S2<FinTT2; i_S2++){
3236 IntPolyh_Triangle &Triangle2 = TTriangles2[i_S2];
3237 if ((Triangle2.IndiceIntersectionPossible() != 0) &&
3238 (Triangle2.GetFleche() >= 0.)) {
3239 //If a triangle is dead or not in BSB, comparison is not possible
3240 Standard_Integer iDeg1, iDeg2, iDeg3, iDeg;
3242 const IntPolyh_Point& P1=TPoints1[Triangle1.FirstPoint()];
3243 const IntPolyh_Point& P2=TPoints1[Triangle1.SecondPoint()];
3244 const IntPolyh_Point& P3=TPoints1[Triangle1.ThirdPoint()];
3245 iDeg1=(P1.Degenerated()) ? 1 : 0;
3246 iDeg2=(P2.Degenerated()) ? 1 : 0;
3247 iDeg3=(P3.Degenerated()) ? 1 : 0;
3248 iDeg=iDeg1+iDeg2+iDeg3;
3253 const IntPolyh_Point& Q1=TPoints2[Triangle2.FirstPoint()];
3254 const IntPolyh_Point& Q2=TPoints2[Triangle2.SecondPoint()];
3255 const IntPolyh_Point& Q3=TPoints2[Triangle2.ThirdPoint()];
3256 iDeg1=(Q1.Degenerated()) ? 1 : 0;
3257 iDeg2=(Q2.Degenerated()) ? 1 : 0;
3258 iDeg3=(Q3.Degenerated()) ? 1 : 0;
3259 iDeg=iDeg1+iDeg2+iDeg3;
3264 if (TriContact(P1, P2, P3, Q1, Q2, Q3, CoupleAngle)) {
3265 if (CpteurTab >= NbTTC)
3267 TTrianglesContacts.SetNbItems(CpteurTab);
3270 TTrianglesContacts[CpteurTab].SetCoupleValue(i_S1, i_S2);
3271 TTrianglesContacts[CpteurTab].SetAngleValue(CoupleAngle);
3272 //test TTrianglesContacts[CpteurTab].Dump(CpteurTab);
3274 Triangle1.SetIndiceIntersection(1);//The triangle is cut by another
3275 Triangle2.SetIndiceIntersection(1);
3281 TTrianglesContacts.SetNbItems(CpteurTab);
3286 //=======================================================================
3287 //function : StartPointsCalcul
3288 //purpose : From the array of couples compute all the start
3289 // points and display them on the screen
3290 //=======================================================================
3291 void IntPolyh_MaillageAffinage::StartPointsCalcul() const
3293 const Standard_Integer FinTTC = TTrianglesContacts.NbItems();
3294 // printf("StartPointsCalcul() from IntPolyh_MaillageAffinage.cxx : StartPoints:\n");
3295 for(Standard_Integer ii=0; ii<FinTTC; ii++) {
3296 IntPolyh_StartPoint SP1,SP2;
3297 Standard_Integer T1,T2;
3298 T1=TTrianglesContacts[ii].FirstValue();
3299 T2=TTrianglesContacts[ii].SecondValue();
3300 StartingPointsResearch(T1,T2,SP1,SP2);
3301 if ( (SP1.E1()!=-1)&&(SP1.E2()!=-1) ) SP1.Dump(ii);
3302 if ( (SP2.E1()!=-1)&&(SP2.E2()!=-1) ) SP2.Dump(ii);
3305 //=======================================================================
3306 //function : CheckCoupleAndGetAngle
3308 //=======================================================================
3309 Standard_Integer CheckCoupleAndGetAngle(const Standard_Integer T1,
3310 const Standard_Integer T2,
3311 Standard_Real& Angle,
3312 IntPolyh_ArrayOfCouples &TTrianglesContacts)
3314 Standard_Integer Test=0;
3315 const Standard_Integer FinTTC = TTrianglesContacts.NbItems();
3316 for (Standard_Integer oioi=0; oioi<FinTTC; oioi++) {
3317 IntPolyh_Couple TestCouple = TTrianglesContacts[oioi];
3318 if ( (TestCouple.FirstValue()==T1)&&(TestCouple.AnalyseFlagValue()!=1) ) {
3319 if (TestCouple.SecondValue()==T2) {
3321 TTrianglesContacts[oioi].SetAnalyseFlag(1);
3322 Angle=TTrianglesContacts[oioi].AngleValue();
3329 //=======================================================================
3330 //function : CheckCoupleAndGetAngle2
3332 //=======================================================================
3333 Standard_Integer CheckCoupleAndGetAngle2(const Standard_Integer T1,
3334 const Standard_Integer T2,
3335 const Standard_Integer T11,
3336 const Standard_Integer T22,
3337 Standard_Integer &CT11,
3338 Standard_Integer &CT22,
3339 Standard_Real & Angle,
3340 IntPolyh_ArrayOfCouples &TTrianglesContacts)
3342 ///couple T1 T2 is found in the list
3343 ///T11 and T22 are two other triangles implied in the contact edge edge
3344 /// CT11 couple( T1,T22) and CT22 couple (T2,T11)
3345 /// these couples will be marked if there is a start point
3346 Standard_Integer Test1=0;
3347 Standard_Integer Test2=0;
3348 Standard_Integer Test3=0;
3349 const Standard_Integer FinTTC = TTrianglesContacts.NbItems();
3350 for (Standard_Integer oioi=0; oioi<FinTTC; oioi++) {
3351 IntPolyh_Couple TestCouple = TTrianglesContacts[oioi];
3352 if( (Test1==0)||(Test2==0)||(Test3==0) ) {
3353 if ( (TestCouple.FirstValue()==T1)&&(TestCouple.AnalyseFlagValue()!=1) ) {
3354 if (TestCouple.SecondValue()==T2) {
3356 TTrianglesContacts[oioi].SetAnalyseFlag(1);
3357 Angle=TTrianglesContacts[oioi].AngleValue();
3359 else if (TestCouple.SecondValue()==T22) {
3362 Angle=TTrianglesContacts[oioi].AngleValue();
3365 else if( (TestCouple.FirstValue()==T11)&&(TestCouple.AnalyseFlagValue()!=1) ) {
3366 if (TestCouple.SecondValue()==T2) {
3369 Angle=TTrianglesContacts[oioi].AngleValue();
3378 //=======================================================================
3379 //function : CheckNextStartPoint
3380 //purpose : it is checked if the point is not a top
3381 // then it is stored in one or several valid arrays with
3382 // the proper list number
3383 //=======================================================================
3384 Standard_Integer CheckNextStartPoint(IntPolyh_SectionLine & SectionLine,
3385 IntPolyh_ArrayOfTangentZones & TTangentZones,
3386 IntPolyh_StartPoint & SP,
3387 const Standard_Boolean Prepend)//=Standard_False)
3389 Standard_Integer Test=1;
3390 if( (SP.E1()==-1)||(SP.E2()==-1) ) {
3391 //The tops of triangle are analyzed
3392 //It is checked if they are not in the array TTangentZones
3393 Standard_Integer FinTTZ=TTangentZones.NbItems();
3394 for(Standard_Integer uiui=0; uiui<FinTTZ; uiui++) {
3395 IntPolyh_StartPoint TestSP=TTangentZones[uiui];
3396 if ( (Abs(SP.U1()-TestSP.U1())<MyConfusionPrecision)
3397 &&(Abs(SP.V1()-TestSP.V1())<MyConfusionPrecision) ) {
3398 if ( (Abs(SP.U2()-TestSP.U2())<MyConfusionPrecision)
3399 &&(Abs(SP.V2()-TestSP.V2())<MyConfusionPrecision) ) {
3400 Test=0;//SP is already in the list of tops
3405 if (Test) {//the top does not belong to the list of TangentZones
3406 SP.SetChainList(-1);
3407 TTangentZones[FinTTZ]=SP;
3408 TTangentZones.IncrementNbItems();
3409 Test=0;//the examined point is a top
3414 SectionLine.Prepend(SP);
3416 SectionLine[SectionLine.NbStartPoints()]=SP;
3417 SectionLine.IncrementNbStartPoints();
3421 //if the point is not a top Test=1
3422 //The chain is continued
3425 //=======================================================================
3426 //function : StartPointsChain
3427 //purpose : Loop on the array of couples. Compute StartPoints.
3428 // Try to chain the StartPoints into SectionLines or
3429 // put the point in the ArrayOfTangentZones if
3430 // chaining it, is not possible.
3431 //=======================================================================
3432 Standard_Integer IntPolyh_MaillageAffinage::StartPointsChain
3433 (IntPolyh_ArrayOfSectionLines& TSectionLines,
3434 IntPolyh_ArrayOfTangentZones& TTangentZones)
3436 //Loop on the array of couples filled in the function COMPARE()
3437 const Standard_Integer FinTTC = TTrianglesContacts.NbItems();
3439 //Array of tops of triangles
3440 for(Standard_Integer IndexA=0; IndexA<FinTTC; IndexA++) {
3441 //First couple of triangles.
3442 //It is checked if the couple of triangles has not been already examined.
3443 if(TTrianglesContacts[IndexA].AnalyseFlagValue()!=1) {
3445 Standard_Integer SectionLineIndex=TSectionLines.NbItems();
3446 // fill last section line if still empty (eap)
3447 if (SectionLineIndex > 0
3449 TSectionLines[SectionLineIndex-1].NbStartPoints() == 0)
3450 SectionLineIndex -= 1;
3452 TSectionLines.IncrementNbItems();
3454 IntPolyh_SectionLine & MySectionLine=TSectionLines[SectionLineIndex];
3455 if (MySectionLine.GetN() == 0) // eap
3456 MySectionLine.Init(10000);//Initialisation of array of StartPoint
3458 Standard_Integer NbPoints=-1;
3459 Standard_Integer T1I, T2I;
3460 T1I = TTrianglesContacts[IndexA].FirstValue();
3461 T2I = TTrianglesContacts[IndexA].SecondValue();
3463 // Start points for the current couple are found
3464 IntPolyh_StartPoint SP1, SP2;
3465 NbPoints=StartingPointsResearch2(T1I,T2I,SP1, SP2);//first calculation
3466 TTrianglesContacts[IndexA].SetAnalyseFlag(1);//the couple is marked
3468 if(NbPoints==1) {// particular case top/triangle or edge/edge
3469 //the start point is input in the array
3470 SP1.SetChainList(SectionLineIndex);
3471 SP1.SetAngle(TTrianglesContacts[IndexA].AngleValue());
3472 //it is checked if the point is not atop of the triangle
3473 if(CheckNextStartPoint(MySectionLine,TTangentZones,SP1)) {
3474 IntPolyh_StartPoint SPNext1;
3475 Standard_Integer TestSP1=0;
3478 IntPolyh_StartPoint SP11;//=SP1;
3479 if(SP1.E1()>=0) { //&&(SP1.E2()!=-1) already tested if the point is not a top
3480 Standard_Integer NextTriangle1=0;
3481 if (TEdges1[SP1.E1()].FirstTriangle()!=T1I) NextTriangle1=TEdges1[SP1.E1()].FirstTriangle();
3482 else NextTriangle1=TEdges1[SP1.E1()].SecondTriangle();
3484 Standard_Real Angle=-2.0;
3485 if (CheckCoupleAndGetAngle(NextTriangle1,T2I,Angle,TTrianglesContacts)) {
3486 //it is checked if the couple exists and is marked
3487 Standard_Integer NbPoints11=0;
3488 NbPoints11=NextStartingPointsResearch2(NextTriangle1,T2I,SP1,SP11);
3489 if (NbPoints11==1) {
3490 SP11.SetChainList(SectionLineIndex);
3491 SP11.SetAngle(Angle);
3493 if(CheckNextStartPoint(MySectionLine,TTangentZones,SP11)) {
3494 Standard_Integer EndChainList=1;
3495 while (EndChainList!=0) {
3496 TestSP1=GetNextChainStartPoint(SP11,SPNext1,MySectionLine,TTangentZones);
3498 SPNext1.SetChainList(SectionLineIndex);
3499 if(CheckNextStartPoint(MySectionLine,TTangentZones,SPNext1))
3501 else EndChainList=0;
3503 else EndChainList=0; //There is no next point
3509 if(NbPoints11>1) {//The point is input in the array TTangentZones
3510 TTangentZones[TTangentZones.NbItems()]=SP11;//default list number = -1
3511 TTangentZones.IncrementNbItems();
3519 else if (SP1.E2()<0){
3522 //chain of the other side
3523 IntPolyh_StartPoint SP12;//=SP1;
3524 if (SP1.E2()>=0) { //&&(SP1.E1()!=-1) already tested
3525 Standard_Integer NextTriangle2;
3526 if (TEdges2[SP1.E2()].FirstTriangle()!=T2I) NextTriangle2=TEdges2[SP1.E2()].FirstTriangle();
3527 else NextTriangle2=TEdges2[SP1.E2()].SecondTriangle();
3529 Standard_Real Angle=-2.0;
3530 if(CheckCoupleAndGetAngle(T1I,NextTriangle2,Angle,TTrianglesContacts)) {
3531 Standard_Integer NbPoints12=0;
3532 NbPoints12=NextStartingPointsResearch2(T1I,NextTriangle2,SP1, SP12);
3533 if (NbPoints12==1) {
3535 SP12.SetChainList(SectionLineIndex);
3536 SP12.SetAngle(Angle);
3537 Standard_Boolean Prepend = Standard_True; // eap
3539 if(CheckNextStartPoint(MySectionLine,TTangentZones,SP12, Prepend)) {
3540 Standard_Integer EndChainList=1;
3541 while (EndChainList!=0) {
3542 TestSP1=GetNextChainStartPoint(SP12,SPNext1,
3543 MySectionLine,TTangentZones,
3546 SPNext1.SetChainList(SectionLineIndex);
3547 if(CheckNextStartPoint(MySectionLine,TTangentZones,SPNext1,Prepend))
3549 else EndChainList=0;
3551 else EndChainList=0; //there is no next point
3556 if(NbPoints12>1) {//The points are input in the array TTangentZones
3557 TTangentZones[TTangentZones.NbItems()]=SP12;//default list number = -1
3558 TTangentZones.IncrementNbItems();
3567 else if(SP1.E1()<0){
3572 else if(NbPoints==2) {
3573 //the start points are input in the array
3574 IntPolyh_StartPoint SPNext2;
3575 Standard_Integer TestSP2=0;
3576 Standard_Integer EndChainList=1;
3578 SP1.SetChainList(SectionLineIndex);
3579 SP1.SetAngle(TTrianglesContacts[IndexA].AngleValue());
3580 if(CheckNextStartPoint(MySectionLine,TTangentZones,SP1)) {
3583 while (EndChainList!=0) {
3584 TestSP2=GetNextChainStartPoint(SP1,SPNext2,MySectionLine,TTangentZones);
3586 SPNext2.SetChainList(SectionLineIndex);
3587 if(CheckNextStartPoint(MySectionLine,TTangentZones,SPNext2))
3589 else EndChainList=0;
3591 else EndChainList=0; //there is no next point
3595 SP2.SetChainList(SectionLineIndex);
3596 SP2.SetAngle(TTrianglesContacts[IndexA].AngleValue());
3597 Standard_Boolean Prepend = Standard_True; // eap
3599 if(CheckNextStartPoint(MySectionLine,TTangentZones,SP2,Prepend)) {
3601 //chain of the other side
3603 while (EndChainList!=0) {
3604 TestSP2=GetNextChainStartPoint(SP2,SPNext2,
3605 MySectionLine,TTangentZones,
3608 SPNext2.SetChainList(SectionLineIndex);
3609 if(CheckNextStartPoint(MySectionLine,TTangentZones,SPNext2,Prepend))
3611 else EndChainList=0;
3613 else EndChainList=0; //there is no next point
3618 else if( (NbPoints>2)&&(NbPoints<7) ) {
3619 //More than two start points
3620 //the start point is input in the table
3621 SP1.SetChainList(SectionLineIndex);
3622 CheckNextStartPoint(MySectionLine,TTangentZones,SP1);
3633 //=======================================================================
3634 //function : GetNextChainStartPoint
3635 //purpose : Mainly used by StartPointsChain(), this function
3636 // try to compute the next StartPoint.
3637 // GetNextChainStartPoint is used only if it is known that there are 2 contact points
3638 //=======================================================================
3639 Standard_Integer IntPolyh_MaillageAffinage::GetNextChainStartPoint
3640 (const IntPolyh_StartPoint & SP,
3641 IntPolyh_StartPoint & SPNext,
3642 IntPolyh_SectionLine & MySectionLine,
3643 IntPolyh_ArrayOfTangentZones & TTangentZones,
3644 const Standard_Boolean Prepend)
3646 Standard_Integer NbPoints=0;
3647 if( (SP.E1()>=0)&&(SP.E2()==-2) ) {
3648 //case if the point is on edge of T1
3649 Standard_Integer NextTriangle1;
3650 if (TEdges1[SP.E1()].FirstTriangle()!=SP.T1()) NextTriangle1=TEdges1[SP.E1()].FirstTriangle();
3652 NextTriangle1=TEdges1[SP.E1()].SecondTriangle();
3653 //If is checked if two triangles intersect
3654 Standard_Real Angle= -2.0;
3655 if (CheckCoupleAndGetAngle(NextTriangle1,SP.T2(),Angle,TTrianglesContacts)) {
3656 NbPoints=NextStartingPointsResearch2(NextTriangle1,SP.T2(),SP,SPNext);
3659 CheckNextStartPoint(MySectionLine,TTangentZones,SPNext,Prepend);
3666 SPNext.SetAngle(Angle);
3668 else NbPoints=0;//this couple does not intersect
3670 else if( (SP.E1()==-2)&&(SP.E2()>=0) ) {
3671 //case if the point is on edge of T2
3672 Standard_Integer NextTriangle2;
3673 if (TEdges2[SP.E2()].FirstTriangle()!=SP.T2()) NextTriangle2=TEdges2[SP.E2()].FirstTriangle();
3675 NextTriangle2=TEdges2[SP.E2()].SecondTriangle();
3676 Standard_Real Angle= -2.0;
3677 if (CheckCoupleAndGetAngle(SP.T1(),NextTriangle2,Angle,TTrianglesContacts)) {
3678 NbPoints=NextStartingPointsResearch2(SP.T1(),NextTriangle2,SP,SPNext);
3681 CheckNextStartPoint(MySectionLine,TTangentZones,SPNext,Prepend);
3688 SPNext.SetAngle(Angle);
3692 else if( (SP.E1()==-2)&&(SP.E2()==-2) ) {
3693 ///no edge is touched or cut
3697 else if( (SP.E1()>=0)&&(SP.E2()>=0) ) {
3698 ///the point is located on two edges
3699 Standard_Integer NextTriangle1;
3700 Standard_Integer CpleT11=-1;
3701 Standard_Integer CpleT22=-1;
3702 if (TEdges1[SP.E1()].FirstTriangle()!=SP.T1()) NextTriangle1=TEdges1[SP.E1()].FirstTriangle();
3704 NextTriangle1=TEdges1[SP.E1()].SecondTriangle();
3705 Standard_Integer NextTriangle2;
3706 if (TEdges2[SP.E2()].FirstTriangle()!=SP.T2()) NextTriangle2=TEdges2[SP.E2()].FirstTriangle();
3708 NextTriangle2=TEdges2[SP.E2()].SecondTriangle();
3709 Standard_Real Angle= -2.0;
3710 if (CheckCoupleAndGetAngle2(NextTriangle1,NextTriangle2,
3711 SP.T1(),SP.T2(),CpleT11,CpleT22,
3712 Angle,TTrianglesContacts)) {
3713 NbPoints=NextStartingPointsResearch2(NextTriangle1,NextTriangle2,SP,SPNext);
3716 ///The new point is checked
3717 if(CheckNextStartPoint(MySectionLine,TTangentZones,SPNext,Prepend)>0) {
3726 SPNext.SetAngle(Angle);
3727 //The couples (Ti,Tj) (Ti',Tj') are marked
3728 if (CpleT11>=0) TTrianglesContacts[CpleT11].SetAnalyseFlag(1);
3732 if (CpleT22>=0) TTrianglesContacts[CpleT22].SetAnalyseFlag(1);
3740 else if( (SP.E1()==-1)||(SP.E2()==-1) ) {
3741 ///the points are tops of triangle
3742 ///the point is atored in an intermediary array
3746 //=======================================================================
3747 //function : GetArrayOfPoints
3749 //=======================================================================
3750 const IntPolyh_ArrayOfPoints& IntPolyh_MaillageAffinage::GetArrayOfPoints
3751 (const Standard_Integer SurfID)const
3757 //=======================================================================
3758 //function : GetArrayOfEdges
3760 //=======================================================================
3761 const IntPolyh_ArrayOfEdges& IntPolyh_MaillageAffinage::GetArrayOfEdges
3762 (const Standard_Integer SurfID)const
3768 //=======================================================================
3769 //function : GetArrayOfTriangles
3771 //=======================================================================
3772 const IntPolyh_ArrayOfTriangles&
3773 IntPolyh_MaillageAffinage::GetArrayOfTriangles
3774 (const Standard_Integer SurfID)const{
3776 return(TTriangles1);
3777 return(TTriangles2);
3780 //=======================================================================
3783 //=======================================================================
3784 Bnd_Box IntPolyh_MaillageAffinage::GetBox(const Standard_Integer SurfID) const
3791 //=======================================================================
3792 //function : GetBoxDraw
3794 //=======================================================================
3795 void IntPolyh_MaillageAffinage::GetBoxDraw(const Standard_Integer SurfID)const
3797 Standard_Real x0,y0,z0,x1,y1,z1;
3799 MyBox1.Get(x0,y0,z0,x1,y1,z1);
3802 MyBox2.Get(x0,y0,z0,x1,y1,z1);
3805 //=======================================================================
3806 //function : GetArrayOfCouples
3808 //=======================================================================
3809 IntPolyh_ArrayOfCouples &IntPolyh_MaillageAffinage::GetArrayOfCouples()
3811 return TTrianglesContacts;
3813 //=======================================================================
3814 //function : SetEnlargeZone
3816 //=======================================================================
3817 void IntPolyh_MaillageAffinage::SetEnlargeZone(Standard_Boolean& EnlargeZone)
3819 myEnlargeZone = EnlargeZone;
3821 //=======================================================================
3822 //function : GetEnlargeZone
3824 //=======================================================================
3825 Standard_Boolean IntPolyh_MaillageAffinage::GetEnlargeZone() const
3827 return myEnlargeZone;
3829 //=======================================================================
3830 //function : DegeneratedIndex
3832 //=======================================================================
3833 void DegeneratedIndex(const TColStd_Array1OfReal& aXpars,
3834 const Standard_Integer aNbX,
3835 const Handle(Adaptor3d_HSurface)& aS,
3836 const Standard_Integer aIsoDirection,
3837 Standard_Integer& aI1,
3838 Standard_Integer& aI2)
3841 Standard_Boolean bDegX1, bDegX2;
3842 Standard_Real aDegX1, aDegX2, aTol2, aX;
3846 aTol2=MyTolerance*MyTolerance;
3848 if (aIsoDirection==1){ // V=const
3849 bDegX1=IsDegenerated(aS, 1, aTol2, aDegX1);
3850 bDegX2=IsDegenerated(aS, 2, aTol2, aDegX2);
3852 else if (aIsoDirection==2){ // U=const
3853 bDegX1=IsDegenerated(aS, 3, aTol2, aDegX1);
3854 bDegX2=IsDegenerated(aS, 4, aTol2, aDegX2);
3860 if (!(bDegX1 || bDegX2)) {
3864 for(i=1; i<=aNbX; ++i) {
3867 if (fabs(aX-aDegX1) < MyTolerance) {
3872 if (fabs(aX-aDegX2) < MyTolerance) {
3878 //=======================================================================
3879 //function : IsDegenerated
3881 //=======================================================================
3882 Standard_Boolean IsDegenerated(const Handle(Adaptor3d_HSurface)& aS,
3883 const Standard_Integer aIndex,
3884 const Standard_Real aTol2,
3885 Standard_Real& aDegX)
3887 Standard_Boolean bRet;
3888 Standard_Integer i, aNbP;
3889 Standard_Real aU, dU, aU1, aU2, aV, dV, aV1, aV2, aD2;
3892 bRet=Standard_False;
3896 aU1=aS->FirstUParameter();
3897 aU2=aS->LastUParameter();
3898 aV1=aS->FirstVParameter();
3899 aV2=aS->LastVParameter();
3901 if (aIndex<3) { // V=const
3906 dU=(aU2-aU1)/(aNbP-1);
3908 aP1=aS->Value(aU, aV);
3909 for (i=1; i<aNbP; ++i) {
3914 aP2=aS->Value(aU, aV);
3915 aD2=aP1.SquareDistance(aP2);
3929 dV=(aV2-aV1)/(aNbP-1);
3931 aP1=aS->Value(aU, aV);
3932 for (i=1; i<aNbP; ++i) {
3937 aP2=aS->Value(aU, aV);
3938 aD2=aP1.SquareDistance(aP2);
3950 //=======================================================================
3951 //function : EnlargeZone
3953 //=======================================================================
3954 void EnlargeZone(const Handle(Adaptor3d_HSurface)& MaSurface,
3960 if(MaSurface->GetType() == GeomAbs_BSplineSurface ||
3961 MaSurface->GetType() == GeomAbs_BezierSurface) {
3962 if((!MaSurface->IsUClosed() && !MaSurface->IsUPeriodic()) &&
3963 (Abs(u0) < 1.e+100 && Abs(u1) < 1.e+100) ) {
3964 Standard_Real delta_u = 0.01*Abs(u1 - u0);
3968 if((!MaSurface->IsVClosed() && !MaSurface->IsVPeriodic()) &&
3969 (Abs(v0) < 1.e+100 && Abs(v1) < 1.e+100) ) {
3970 Standard_Real delta_v = 0.01*Abs(v1 - v0);
3978 #include <TopoDS_Shape.hxx>
3979 #include <Poly_Triangulation.hxx>
3980 #include <TColgp_Array1OfPnt.hxx>
3981 #include <Poly_Array1OfTriangle.hxx>
3982 #include <BRep_TFace.hxx>
3983 #include <TopoDS_Face.hxx>
3985 //=======================================================================
3986 //function : TriangleShape
3987 //purpose : shape with triangulation containing triangles
3988 //=======================================================================
3989 static TopoDS_Shape TriangleShape(const IntPolyh_ArrayOfTriangles & TTriangles,
3990 const IntPolyh_ArrayOfPoints & TPoints)
3993 if (TPoints.NbItems() < 1 || TTriangles.NbItems() < 1) return aFace;
3995 Handle(Poly_Triangulation) aPTriangulation =
3996 new Poly_Triangulation(TPoints.NbItems(),TTriangles.NbItems(),Standard_False);
3997 TColgp_Array1OfPnt & aPNodes = aPTriangulation->ChangeNodes();
3998 Poly_Array1OfTriangle & aPTrialgles = aPTriangulation->ChangeTriangles();
4001 for (i=0; i<TPoints.NbItems(); i++) {
4002 const IntPolyh_Point& P = TPoints[i];
4003 aPNodes(i+1).SetCoord(P.X(), P.Y(), P.Z());
4005 for (i=0; i<TTriangles.NbItems(); i++) {
4006 const IntPolyh_Triangle& T = TTriangles[i];
4007 aPTrialgles(i+1).Set(T.FirstPoint()+1, T.SecondPoint()+1, T.ThirdPoint()+1);
4010 Handle(BRep_TFace) aTFace = new BRep_TFace;
4011 aTFace->Triangulation(aPTriangulation);
4012 aFace.TShape(aTFace);
4017 //#define MyTolerance 10.0e-7
4018 //#define MyConfusionPrecision 10.0e-12
4019 //#define SquareMyConfusionPrecision 10.0e-24