1 // Created on: 1993-02-03
2 // Created by: Laurent BUCHARD
3 // Copyright (c) 1993-1999 Matra Datavision
4 // Copyright (c) 1999-2014 OPEN CASCADE SAS
6 // This file is part of Open CASCADE Technology software library.
8 // This library is free software; you can redistribute it and/or modify it under
9 // the terms of the GNU Lesser General Public License version 2.1 as published
10 // by the Free Software Foundation, with special exception defined in the file
11 // OCCT_LGPL_EXCEPTION.txt. Consult the file LICENSE_LGPL_21.txt included in OCCT
12 // distribution for complete text of the license and disclaimer of any warranty.
14 // Alternatively, this file may be used under the terms of Open CASCADE
15 // commercial license or contractual agreement.
18 #include <Adaptor3d_HSurface.hxx>
19 #include <Bnd_Array1OfBox.hxx>
20 #include <Bnd_Box.hxx>
24 #include <IntPatch_HInterTool.hxx>
25 #include <IntPatch_Polyhedron.hxx>
26 #include <Precision.hxx>
27 #include <Standard_ConstructionError.hxx>
28 #include <Standard_OutOfRange.hxx>
29 #include <TColgp_Array2OfPnt.hxx>
30 #include <TColStd_Array2OfReal.hxx>
35 #define LONGUEUR_MINI_EDGE_TRIANGLE 1e-14
36 #define DEFLECTION_COEFF 1.1
39 //================================================================================
40 static Standard_Integer NbPOnU (const Handle(Adaptor3d_HSurface)& S)
42 const Standard_Real u0 = S->FirstUParameter();
43 const Standard_Real u1 = S->LastUParameter();
44 const Standard_Integer nbpu = IntPatch_HInterTool::NbSamplesU(S,u0,u1);
45 return (nbpu>NBMAXUV? NBMAXUV : nbpu);
47 //================================================================================
48 static Standard_Integer NbPOnV (const Handle(Adaptor3d_HSurface)& S)
50 const Standard_Real v0 = S->FirstVParameter();
51 const Standard_Real v1 = S->LastVParameter();
52 const Standard_Integer nbpv = IntPatch_HInterTool::NbSamplesV(S,v0,v1);
53 return (nbpv>NBMAXUV? NBMAXUV : nbpv);
56 //=======================================================================
59 //=======================================================================
60 void IntPatch_Polyhedron::Destroy()
62 gp_Pnt *CMyPnts = (gp_Pnt *)C_MyPnts; if(C_MyPnts) delete [] CMyPnts;
63 Standard_Real *CMyU = (Standard_Real *)C_MyU; if(C_MyU) delete [] CMyU;
64 Standard_Real *CMyV = (Standard_Real *)C_MyV; if(C_MyV) delete [] CMyV;
65 C_MyPnts=C_MyU=C_MyV=NULL;
68 //=======================================================================
69 //function : IntPatch_Polyhedron
71 //=======================================================================
72 IntPatch_Polyhedron::IntPatch_Polyhedron (const Handle(Adaptor3d_HSurface)& Surface)
73 : TheDeflection(Epsilon(100.)),
74 nbdeltaU(NbPOnU(Surface)),
75 nbdeltaV(NbPOnV(Surface)),
76 C_MyPnts(NULL),C_MyU(NULL),C_MyV(NULL),
77 UMinSingular(IntPatch_HInterTool::SingularOnVMin(Surface)),
78 UMaxSingular(IntPatch_HInterTool::SingularOnVMin(Surface)),
79 VMinSingular(IntPatch_HInterTool::SingularOnVMin(Surface)),
80 VMaxSingular(IntPatch_HInterTool::SingularOnVMin(Surface))
82 const Standard_Integer t = (nbdeltaU+1)*(nbdeltaV+1)+1;
83 gp_Pnt *CMyPnts = new gp_Pnt[t];
84 Standard_Real *CMyU = new Standard_Real[t];
85 Standard_Real *CMyV = new Standard_Real[t];
90 const Standard_Real u0 = Surface->FirstUParameter();
91 const Standard_Real u1 = Surface->LastUParameter();
92 const Standard_Real v0 = Surface->FirstVParameter();
93 const Standard_Real v1 = Surface->LastVParameter();
95 const Standard_Real U1mU0sNbdeltaU = (u1-u0)/(Standard_Real)nbdeltaU;
96 const Standard_Real V1mV0sNbdeltaV = (v1-v0)/(Standard_Real)nbdeltaV;
100 Standard_Integer i1, i2, Index=1;
101 for (i1 = 0, U = u0; i1 <= nbdeltaU; i1++, U+= U1mU0sNbdeltaU) {
102 for (i2 = 0, V = v0; i2 <= nbdeltaV; i2++, V+= V1mV0sNbdeltaV ) {
112 Standard_Real tol=0.0;
113 const Standard_Integer nbtriangles = NbTriangles();
114 for (i1=1; i1<=nbtriangles; i1++) {
115 const Standard_Real tol1 = DeflectionOnTriangle(Surface,i1);
116 if(tol1>tol) tol=tol1;
119 tol*=DEFLECTION_COEFF;
121 DeflectionOverEstimation(tol);
125 //=======================================================================
126 //function : IntPatch_Polyhedron
128 //=======================================================================
129 IntPatch_Polyhedron::IntPatch_Polyhedron (const Handle(Adaptor3d_HSurface)& Surface,
130 const Standard_Integer nbu,
131 const Standard_Integer nbv)
132 : TheDeflection(Epsilon(100.)),
135 C_MyPnts(NULL),C_MyU(NULL),C_MyV(NULL),
136 UMinSingular(IntPatch_HInterTool::SingularOnVMin(Surface)),
137 UMaxSingular(IntPatch_HInterTool::SingularOnVMin(Surface)),
138 VMinSingular(IntPatch_HInterTool::SingularOnVMin(Surface)),
139 VMaxSingular(IntPatch_HInterTool::SingularOnVMin(Surface))
141 const Standard_Integer t = (nbdeltaU+1)*(nbdeltaV+1)+1;
142 gp_Pnt *CMyPnts = new gp_Pnt[t];
143 Standard_Real *CMyU = new Standard_Real[t];
144 Standard_Real *CMyV = new Standard_Real[t];
149 const Standard_Real u0 = Surface->FirstUParameter();
150 const Standard_Real u1 = Surface->LastUParameter();
151 const Standard_Real v0 = Surface->FirstVParameter();
152 const Standard_Real v1 = Surface->LastVParameter();
154 const Standard_Real U1mU0sNbdeltaU = (u1-u0)/(Standard_Real)nbdeltaU;
155 const Standard_Real V1mV0sNbdeltaV = (v1-v0)/(Standard_Real)nbdeltaV;
159 Standard_Integer i1, i2, Index=1;
160 for (i1 = 0, U = u0; i1 <= nbdeltaU; i1++, U+= U1mU0sNbdeltaU) {
161 for (i2 = 0, V = v0; i2 <= nbdeltaV; i2++, V+= V1mV0sNbdeltaV ) {
171 Standard_Real tol=0.0;
172 const Standard_Integer nbtriangles = NbTriangles();
173 for (i1=1; i1<=nbtriangles; i1++) {
174 const Standard_Real tol1 = DeflectionOnTriangle(Surface,i1);
175 if(tol1>tol) tol=tol1;
178 tol*=DEFLECTION_COEFF;
180 DeflectionOverEstimation(tol);
185 //=======================================================================
186 //function : DeflectionOnTriangle
188 //=======================================================================
190 Standard_Real IntPatch_Polyhedron::DeflectionOnTriangle
191 (const Handle(Adaptor3d_HSurface)& Surface,
192 const Standard_Integer Triang) const
194 Standard_Integer i1,i2,i3;
196 Triangle(Triang,i1,i2,i3);
197 //-- Calcul de l eqution du plan
198 Standard_Real u1,v1,u2,v2,u3,v3;
200 P1 = Point(i1,u1,v1);
201 P2 = Point(i2,u2,v2);
202 P3 = Point(i3,u3,v3);
203 if(P1.SquareDistance(P2)<=LONGUEUR_MINI_EDGE_TRIANGLE) return(0);
204 if(P1.SquareDistance(P3)<=LONGUEUR_MINI_EDGE_TRIANGLE) return(0);
205 if(P2.SquareDistance(P3)<=LONGUEUR_MINI_EDGE_TRIANGLE) return(0);
206 gp_XYZ XYZ1=P2.XYZ()-P1.XYZ();
207 gp_XYZ XYZ2=P3.XYZ()-P2.XYZ();
208 gp_XYZ XYZ3=P1.XYZ()-P3.XYZ();
209 gp_Vec NormalVector((XYZ1^XYZ2)+(XYZ2^XYZ3)+(XYZ3^XYZ1));
210 Standard_Real aNormLen = NormalVector.Magnitude();
211 if (aNormLen < gp::Resolution()) {
215 NormalVector.Divide(aNormLen);
216 //-- Calcul du point u,v au centre du triangle
217 Standard_Real u = (u1+u2+u3)/3.0;
218 Standard_Real v = (v1+v2+v3)/3.0;
219 gp_Vec P1P(P1,Surface->Value(u,v));
220 return(Abs(P1P.Dot(NormalVector)));
223 //=======================================================================
224 //function : Parameters
226 //=======================================================================
227 void IntPatch_Polyhedron::Parameters( const Standard_Integer Index
229 ,Standard_Real &V) const
231 U = ((Standard_Real *)C_MyU)[Index];
232 V = ((Standard_Real *)C_MyV)[Index];
235 //=======================================================================
236 //function : DeflectionOverEstimation
238 //=======================================================================
239 void IntPatch_Polyhedron::DeflectionOverEstimation(const Standard_Real flec)
242 TheDeflection=0.0001;
243 TheBnd.Enlarge(0.0001);
247 TheBnd.Enlarge(flec);
251 //=======================================================================
252 //function : DeflectionOverEstimation
254 //=======================================================================
255 Standard_Real IntPatch_Polyhedron::DeflectionOverEstimation() const
257 return TheDeflection;
260 //=======================================================================
261 //function : Bounding
263 //=======================================================================
264 const Bnd_Box& IntPatch_Polyhedron::Bounding() const
269 //=======================================================================
270 //function : FillBounding
272 //=======================================================================
273 void IntPatch_Polyhedron::FillBounding()
275 TheComponentsBnd=new Bnd_HArray1OfBox(1, NbTriangles());
277 Standard_Integer p1, p2, p3;
278 Standard_Integer nbtriangles = NbTriangles();
279 for (Standard_Integer iTri=1; iTri<=nbtriangles; iTri++) {
280 Triangle(iTri, p1, p2, p3);
282 const gp_Pnt& P1 = Point(p1);
283 const gp_Pnt& P2 = Point(p2);
284 const gp_Pnt& P3 = Point(p3);
285 if(P1.SquareDistance(P2)>LONGUEUR_MINI_EDGE_TRIANGLE) {
286 if(P1.SquareDistance(P3)>LONGUEUR_MINI_EDGE_TRIANGLE) {
287 if(P2.SquareDistance(P3)>LONGUEUR_MINI_EDGE_TRIANGLE) {
294 Boite.Enlarge(TheDeflection);
295 TheComponentsBnd->SetValue(iTri,Boite);
299 //=======================================================================
300 //function : ComponentsBounding
302 //=======================================================================
303 const Handle(Bnd_HArray1OfBox)& IntPatch_Polyhedron::ComponentsBounding () const
305 return TheComponentsBnd;
308 //=======================================================================
309 //function : NbTriangles
311 //=======================================================================
312 Standard_Integer IntPatch_Polyhedron::NbTriangles () const
314 return nbdeltaU*nbdeltaV*2;
317 //=======================================================================
318 //function : NbPoints
320 //=======================================================================
321 Standard_Integer IntPatch_Polyhedron::NbPoints () const
323 return (nbdeltaU+1)*(nbdeltaV+1);
326 //=======================================================================
327 //function : TriConnex
329 //=======================================================================
330 Standard_Integer IntPatch_Polyhedron::TriConnex (const Standard_Integer Triang,
331 const Standard_Integer Pivot,
332 const Standard_Integer Pedge,
333 Standard_Integer& TriCon,
334 Standard_Integer& OtherP) const {
336 Standard_Integer Pivotm1 = Pivot-1;
337 Standard_Integer nbdeltaVp1 = nbdeltaV+1;
338 Standard_Integer nbdeltaVm2 = nbdeltaV + nbdeltaV;
340 // Pivot position in the MaTriangle :
341 Standard_Integer ligP = Pivotm1/nbdeltaVp1;
342 Standard_Integer colP = Pivotm1 - ligP * nbdeltaVp1;
344 // Point sur Edge position in the MaTriangle and edge typ :
345 Standard_Integer ligE = 0, colE = 0, typE = 0;
347 ligE= (Pedge-1)/nbdeltaVp1;
348 colE= (Pedge-1) - (ligE * nbdeltaVp1);
350 if (ligP==ligE) typE=1;
352 else if (colP==colE) typE=2;
360 // Triangle position General case :
361 Standard_Integer linT = 0, colT = 0;
362 Standard_Integer linO = 0, colO = 0;
363 Standard_Integer t,tt;
365 t = (Triang-1)/(nbdeltaVm2);
366 tt= (Triang-1)-t*nbdeltaVm2;
376 if (colT==ligP+ligP) {
389 case 1: // Horizontal
393 colO=(colP>colE)? colP : colE; //--colO=Max(colP, colE);
398 colO=(colP<colE)? colP : colE; //--colO=Min(colP, colE);
402 if (colT==(colP+colP)) {
404 linO=(ligP>ligE)? ligP : ligE; //--linO=Max(ligP, ligE);
409 linO=(ligP<ligE)? ligP : ligE; //--linO=Min(ligP, ligE);
416 linO=(ligP>ligE)? ligP : ligE; //--linO=Max(ligP, ligE);
417 colO=(colP<colE)? colP : colE; //--colO=Min(colP, colE);
421 linO=(ligP<ligE)? ligP : ligE; //--linO=Min(ligP, ligE);
422 colO=(colP>colE)? colP : colE; //--colO=Max(colP, colE);
428 // Unknown Triangle position :
431 linT=(1>ligP)? 1 : ligP; //--linT=Max(1, ligP);
432 colT=(1>(colP+colP))? 1 : (colP+colP); //--colT=Max(1, colP+colP);
433 if (ligP==0) linO=ligP+1;
438 // Known edge We take the left or down connectivity :
440 case 1: // Horizontal
442 colT=(colP>colE)? colP : colE; //--colT=Max(colP,colE);
445 colO=(colP>colE)? colP : colE; //--colO=Max(colP,colE);
448 linT=(ligP>ligE)? ligP : ligE; //--linT=Max(ligP, ligE);
450 linO=(ligP<ligE)? ligP : ligE; //--linO=Min(ligP, ligE);
454 linT=(ligP>ligE)? ligP : ligE; //--linT=Max(ligP, ligE);
456 linO=(ligP>ligE)? ligP : ligE; //--linO=Max(ligP, ligE);
457 colO=(colP<colE)? colP : colE; //--colO=Min(colP, colE);
463 TriCon=(linT-1)*nbdeltaVm2 + colT;
468 if (colO<0) {colO=0;linO=1;}
469 else if (colO>nbdeltaV) {colO=nbdeltaV;linO=1;}
472 else if (linT>nbdeltaU) {
475 if (colO<0) {colO=0;linO=nbdeltaU-1;}
476 else if (colO>nbdeltaV) {colO=nbdeltaV;linO=nbdeltaU-1;}
483 if (linO<0) {linO=0;colO=1;}
484 else if (linO>nbdeltaU) {linO=nbdeltaU;colO=1;}
487 else if (colT>nbdeltaV) {
490 if (linO<0) {linO=0;colO=nbdeltaV-1;}
491 else if (linO>nbdeltaU) {linO=nbdeltaU;colO=nbdeltaV-1;}
495 OtherP=linO*nbdeltaVp1 + colO+1;
498 //----------------------------------------------------
499 //-- Detection des cas ou le triangle retourne est
500 //-- invalide. Dans ce cas, on retourne le triangle
501 //-- suivant par un nouvel appel a TriConnex.
503 //-- Si En entree : Point(Pivot)==Point(Pedge)
504 //-- Alors on retourne OtherP a 0
505 //-- et Tricon = Triangle
507 if(Point(Pivot).SquareDistance(Point(Pedge))<=LONGUEUR_MINI_EDGE_TRIANGLE) {
511 cout<<" Probleme ds IntCurveSurface_Polyhedron : Pivot et PEdge Confondus "<<endl;
515 if(Point(OtherP).SquareDistance(Point(Pedge))<=LONGUEUR_MINI_EDGE_TRIANGLE) {
517 cout<<" Probleme ds IntCurveSurface_Polyhedron : OtherP et PEdge Confondus "<<endl;
519 return(0); //-- BUG NON CORRIGE ( a revoir le role de nbdeltaU et nbdeltaV)
520 // Standard_Integer TempTri,TempOtherP;
522 // TempOtherP = OtherP;
523 // return(TriConnex(TempTri,Pivot,TempOtherP,TriCon,OtherP));
530 //=======================================================================
531 //function : PlaneEquation
533 //=======================================================================
535 void IntPatch_Polyhedron::PlaneEquation (const Standard_Integer Triang,
536 gp_XYZ& NormalVector,
537 Standard_Real& PolarDistance) const
539 Standard_Integer i1,i2,i3;
540 Triangle(Triang,i1,i2,i3);
542 gp_XYZ Pointi1(Point(i1).XYZ());
543 gp_XYZ Pointi2(Point(i2).XYZ());
544 gp_XYZ Pointi3(Point(i3).XYZ());
547 gp_XYZ v1= Pointi2 - Pointi1;
548 gp_XYZ v2= Pointi3 - Pointi2;
549 gp_XYZ v3= Pointi1 - Pointi3;
551 if(v1.SquareModulus()<=LONGUEUR_MINI_EDGE_TRIANGLE) { NormalVector.SetCoord(1.0,0.0,0.0); return; }
552 if(v2.SquareModulus()<=LONGUEUR_MINI_EDGE_TRIANGLE) { NormalVector.SetCoord(1.0,0.0,0.0); return; }
553 if(v3.SquareModulus()<=LONGUEUR_MINI_EDGE_TRIANGLE) { NormalVector.SetCoord(1.0,0.0,0.0); return; }
555 NormalVector= (v1^v2)+(v2^v3)+(v3^v1);
556 Standard_Real aNormLen = NormalVector.Modulus();
557 if (aNormLen < gp::Resolution()) {
561 NormalVector.Divide(aNormLen);
562 PolarDistance = NormalVector * Point(i1).XYZ();
565 //=======================================================================
568 //=======================================================================
569 Standard_Boolean IntPatch_Polyhedron::Contain (const Standard_Integer Triang,
570 const gp_Pnt& ThePnt) const
572 Standard_Integer i1,i2,i3;
573 Triangle(Triang,i1,i2,i3);
574 gp_XYZ Pointi1(Point(i1).XYZ());
575 gp_XYZ Pointi2(Point(i2).XYZ());
576 gp_XYZ Pointi3(Point(i3).XYZ());
578 gp_XYZ v1=(Pointi2-Pointi1)^(ThePnt.XYZ()-Pointi1);
579 gp_XYZ v2=(Pointi3-Pointi2)^(ThePnt.XYZ()-Pointi2);
580 gp_XYZ v3=(Pointi1-Pointi3)^(ThePnt.XYZ()-Pointi3);
581 if (v1*v2 >= 0. && v2*v3 >= 0. && v3*v1>=0.)
582 return Standard_True;
584 return Standard_False;
586 //=======================================================================
589 //=======================================================================
591 void IntPatch_Polyhedron::Dump()const
594 //=======================================================================
597 //=======================================================================
598 void IntPatch_Polyhedron::Size(Standard_Integer& nbdu,
599 Standard_Integer& nbdv) const
604 //=======================================================================
605 //function : Triangle
607 //=======================================================================
608 void IntPatch_Polyhedron::Triangle (const Standard_Integer Index,
609 Standard_Integer& P1,
610 Standard_Integer& P2,
611 Standard_Integer& P3) const
613 Standard_Integer line=1+((Index-1)/(nbdeltaV*2));
614 Standard_Integer colon=1+((Index-1)%(nbdeltaV*2));
615 Standard_Integer colpnt=(colon+1)/2;
617 // General formula = (line-1)*(nbdeltaV+1)+colpnt
619 // Position of P1 = MesXYZ(line,colpnt);
620 P1= (line-1)*(nbdeltaV+1) + colpnt;
622 // Position of P2= MesXYZ(line+1,colpnt+((colon-1)%2));
623 P2= line*(nbdeltaV+1) + colpnt+((colon-1)%2);
625 // Position of P3= MesXYZ(line+(colon%2),colpnt+1);
626 P3= (line-1+(colon%2))*(nbdeltaV+1) + colpnt + 1;
627 //-- printf("\nTriangle %4d P1:%4d P2:%4d P3:%4d",Index,P1,P2,P3);
629 //=======================================================================
631 //=======================================================================
632 const gp_Pnt& IntPatch_Polyhedron::Point( const Standard_Integer Index
634 ,Standard_Real& V) const
636 gp_Pnt *CMyPnts = (gp_Pnt *)C_MyPnts;
637 Standard_Real *CMyU = (Standard_Real *)C_MyU;
638 Standard_Real *CMyV = (Standard_Real *)C_MyV;
641 return CMyPnts[Index];
643 //=======================================================================
645 //=======================================================================
646 const gp_Pnt& IntPatch_Polyhedron::Point(const Standard_Integer Index) const {
647 gp_Pnt *CMyPnts = (gp_Pnt *)C_MyPnts;
648 return CMyPnts[Index];
651 //=======================================================================
653 //=======================================================================
654 void IntPatch_Polyhedron::Point (const gp_Pnt& /*p*/,
655 const Standard_Integer /*lig*/,
656 const Standard_Integer /*col*/,
657 const Standard_Real /*u*/,
658 const Standard_Real /*v*/)
660 //printf("\n IntPatch_Polyhedron::Point : Ne dois pas etre appelle\n");
663 //=======================================================================
665 //=======================================================================
666 void IntPatch_Polyhedron::Point (const Standard_Integer Index, gp_Pnt& P) const
668 gp_Pnt *CMyPnts = (gp_Pnt *)C_MyPnts;
671 //=======================================================================