0024530: TKMesh - remove unused package IntPoly
[occt.git] / src / MeshTest / MeshTest.cxx
1 // Created on: 1993-09-22
2 // Created by: Didier PIFFAULT
3 // Copyright (c) 1993-1999 Matra Datavision
4 // Copyright (c) 1999-2014 OPEN CASCADE SAS
5 //
6 // This file is part of Open CASCADE Technology software library.
7 //
8 // This library is free software; you can redistribute it and / or modify it
9 // under the terms of the GNU Lesser General Public version 2.1 as published
10 // by the Free Software Foundation, with special exception defined in the file
11 // OCCT_LGPL_EXCEPTION.txt. Consult the file LICENSE_LGPL_21.txt included in OCCT
12 // distribution for complete text of the license and disclaimer of any warranty.
13 //
14 // Alternatively, this file may be used under the terms of Open CASCADE
15 // commercial license or contractual agreement.
16
17 #include <Standard_Stream.hxx>
18
19 #include <stdio.h>
20
21 #include <MeshTest.ixx>
22
23 #include <MeshTest_DrawableMesh.hxx>
24 #include <TopAbs_ShapeEnum.hxx>
25 #include <TopoDS.hxx>
26 #include <TopoDS_Edge.hxx>
27 #include <TopoDS_Face.hxx>
28 #include <TopoDS_Shape.hxx>
29 #include <TopoDS_Compound.hxx>
30 #include <TopExp_Explorer.hxx>
31 #include <TopTools_ListIteratorOfListOfShape.hxx>
32 #include <DBRep.hxx>
33 #include <BRepTest.hxx>
34 #include <GeometryTest.hxx>
35 #include <BRep_Tool.hxx>
36 #include <BRep_Builder.hxx>
37 #include <Draw_MarkerShape.hxx>
38 #include <Draw_Appli.hxx>
39 #include <Draw.hxx>
40 #include <DrawTrSurf.hxx>
41 #include <BRepMesh_Triangle.hxx>
42 #include <BRepMesh_DataStructureOfDelaun.hxx>
43 #include <BRepMesh_Delaun.hxx>
44 #include <BRepMesh_FastDiscret.hxx>
45 #include <BRepMesh_Array1OfVertexOfDelaun.hxx>
46 #include <BRepMesh_Vertex.hxx>
47 #include <BRepMesh_Edge.hxx>
48 #include <BRepMesh_IncrementalMesh.hxx>
49 #include <TColStd_ListIteratorOfListOfInteger.hxx>
50 #include <TColStd_MapIteratorOfMapOfInteger.hxx>
51 #include <Bnd_Box.hxx>
52 #include <Precision.hxx>
53 #include <Draw_Interpretor.hxx>
54 #include <Geom_Plane.hxx>
55 #include <Geom_Surface.hxx>
56 #include <Draw_Marker3D.hxx>
57 #include <Draw_Segment2D.hxx>
58
59 #include <GCPnts_UniformAbscissa.hxx>
60 #include <GeomAdaptor_Curve.hxx>
61 #include <Geom_Curve.hxx>
62 #include <Extrema_LocateExtPC.hxx>
63
64 #include <TopLoc_Location.hxx>
65 #include <gp_Trsf.hxx>
66 #include <Poly_Triangulation.hxx>
67 #include <Poly_Connect.hxx>
68 #include <TColgp_Array1OfPnt2d.hxx>
69 #include <TColStd_HArray1OfInteger.hxx>
70 #include <TopExp_Explorer.hxx>
71 #include <gp_Pln.hxx>
72
73 #include <PLib.hxx>
74 #include <AppCont_ContMatrices.hxx>
75 #include <math_Vector.hxx>
76 #include <math_Matrix.hxx>
77 #include <math.hxx>
78
79 #include <CSLib_DerivativeStatus.hxx>
80 #include <CSLib.hxx>
81 #include <BRepAdaptor_Surface.hxx>
82 #include <Bnd_Box.hxx>
83 #include <BRepBndLib.hxx>
84
85
86 //epa Memory leaks test
87 #include <BRepBuilderAPI_MakePolygon.hxx>
88 #include <TopoDS_Wire.hxx>
89 #include <BRepBuilderAPI_MakeFace.hxx>
90 #include <BRepTools.hxx>
91
92 //OAN: for triepoints
93 #include <BRepBuilderAPI_MakeVertex.hxx>
94 #include <Poly_PolygonOnTriangulation.hxx>
95 #include <TopTools_MapIteratorOfMapOfShape.hxx>
96
97 #ifdef WNT
98 Standard_IMPORT Draw_Viewer dout;
99 #endif
100
101 #define MAX2(X, Y)      (  Abs(X) > Abs(Y)? Abs(X) : Abs(Y) )
102 #define MAX3(X, Y, Z)   ( MAX2 ( MAX2(X,Y) , Z) )
103
104
105
106 #define ONETHIRD 0.333333333333333333333333333333333333333333333333333333333333
107 #define TWOTHIRD 0.666666666666666666666666666666666666666666666666666666666666
108
109 #ifdef DEB_MESH_CHRONO
110 #include <OSD_Chronometer.hxx>
111 Standard_Integer D0Control, D0Internal, D0Unif, D0Edges, NbControls;
112 OSD_Chronometer chTotal, chInternal, chControl, chUnif, chAddPoint;
113 OSD_Chronometer chEdges, chMaillEdges, chEtuInter, chLastControl, chStock;
114 OSD_Chronometer chAdd11, chAdd12, chAdd2, chUpdate, chPointValid;
115 OSD_Chronometer chIsos, chPointsOnIsos;
116 #endif
117
118 //=======================================================================
119 //function : incrementalmesh
120 //purpose  : 
121 //=======================================================================
122
123 static Standard_Integer incrementalmesh(Draw_Interpretor& di, Standard_Integer nbarg, const char** argv)
124 {
125   if (nbarg < 3) {
126     di << " use incmesh shape deflection [inParallel (0/1) : 0 by default]\n";
127     return 0;
128   }
129
130   TopoDS_Shape aShape = DBRep::Get(argv[1]);
131   if (aShape.IsNull()) {
132     di << " null shapes is not allowed here\n";
133     return 0;
134   }
135   Standard_Real aDeflection = Draw::Atof(argv[2]);
136
137   Standard_Boolean isInParallel = Standard_False;
138   if (nbarg == 4) {
139         isInParallel = Draw::Atoi(argv[3]) == 1;
140   }
141   di << "Incremental Mesh, multi-threading "
142     << (isInParallel ? "ON\n" : "OFF\n");
143   
144   BRepMesh_IncrementalMesh MESH(aShape, aDeflection, Standard_False, 0.5, isInParallel);
145   Standard_Integer statusFlags = MESH.GetStatusFlags();  
146
147   di << "Meshing statuses: ";
148
149   if( !statusFlags )
150   {
151     di << "NoError";
152   }
153   else
154   {
155     Standard_Integer i;
156     for( i = 0; i < 4; i++ )
157     {
158       if( (statusFlags >> i) & (Standard_Integer)1 )
159       {
160         switch(i+1)
161         {
162           case 1:
163             di << "OpenWire ";
164             break;
165           case 2:
166             di << "SelfIntersectingWire ";
167             break;
168           case 3:
169             di << "Failure ";
170             break;
171           case 4:
172             di << "ReMesh ";
173             break;
174         }
175       }
176     }
177   }
178
179   return 0;
180 }
181
182 //=======================================================================
183 //function : MemLeakTest
184 //purpose  : 
185 //=======================================================================
186
187 static Standard_Integer MemLeakTest(Draw_Interpretor&, Standard_Integer /*nbarg*/, const char** /*argv*/)
188 {
189   for(int i=0;i<10000;i++)
190   {
191     BRepBuilderAPI_MakePolygon w(gp_Pnt(0,0,0),gp_Pnt(0,100,0),gp_Pnt(20,100,0),gp_Pnt(20,0,0));
192     w.Close();     
193     TopoDS_Wire wireShape( w.Wire());
194     BRepBuilderAPI_MakeFace faceBuilder(wireShape);          
195     TopoDS_Face f( faceBuilder.Face());
196     BRepMesh_IncrementalMesh im(f,1);
197     BRepTools::Clean(f);      
198   }
199   return 0;
200 }
201
202 //=======================================================================
203 //function : fastdiscret
204 //purpose  : 
205 //=======================================================================
206
207 static Standard_Integer fastdiscret(Draw_Interpretor& di, Standard_Integer nbarg, const char** argv)
208 {
209   if (nbarg < 3) return 1;
210
211   TopoDS_Shape S = DBRep::Get(argv[1]);
212   if (S.IsNull()) return 1;
213
214   const Standard_Real d = Draw::Atof(argv[2]);
215
216   Standard_Boolean WithShare = Standard_True;
217   if (nbarg > 3) WithShare = Draw::Atoi(argv[3]);
218
219   Bnd_Box B;
220   BRepBndLib::Add(S,B);
221   BRepMesh_FastDiscret MESH(d,0.5,B,WithShare,Standard_True,Standard_False,Standard_True);
222
223   //Standard_Integer NbIterations = MESH.NbIterations();
224   //if (nbarg > 4) NbIterations = Draw::Atoi(argv[4]);
225   //MESH.NbIterations() = NbIterations;
226
227   di<<"Starting FastDiscret with :"<<"\n";
228   di<<"  Deflection="<<d<<"\n";
229   di<<"  Angle="<<0.5<<"\n";
230   di<<"  SharedMode="<< (Standard_Integer) WithShare<<"\n";
231   //di<<"  NbIterations="<<NbIterations<<"\n";
232
233   Handle(Poly_Triangulation) T;
234   BRep_Builder aBuilder;
235   TopExp_Explorer ex;
236
237   // Clear existing triangulations
238   for (ex.Init(S, TopAbs_FACE); ex.More(); ex.Next())
239     aBuilder.UpdateFace(TopoDS::Face(ex.Current()),T);
240
241   MESH.Perform(S);
242
243   TopoDS_Compound aCompGood, aCompFailed, aCompViolating;
244
245   TopLoc_Location L;
246   Standard_Integer nbtriangles = 0, nbnodes = 0, nbfailed = 0, nbviolating = 0;
247   Standard_Real maxdef = 0.0;
248   for (ex.Init(S, TopAbs_FACE); ex.More(); ex.Next())
249   {
250     T = BRep_Tool::Triangulation(TopoDS::Face(ex.Current()),L);
251     if (T.IsNull())
252     {
253       nbfailed++;
254       if (aCompFailed.IsNull())
255         aBuilder.MakeCompound(aCompFailed);
256       aBuilder.Add(aCompFailed,ex.Current());
257     }
258     else
259     {
260       nbtriangles += T->NbTriangles();
261       nbnodes += T->NbNodes();
262       if (T->Deflection() > maxdef) maxdef = T->Deflection();
263       if (T->Deflection() > d)
264       {
265         nbviolating++;
266         if (aCompViolating.IsNull())
267           aBuilder.MakeCompound(aCompViolating);
268         aBuilder.Add(aCompViolating,ex.Current());
269       }
270       else
271       {
272         if (aCompGood.IsNull())
273           aBuilder.MakeCompound(aCompGood);
274         aBuilder.Add(aCompGood,ex.Current());
275       }
276     }
277   }
278
279   if (!aCompGood.IsNull())
280   {
281     char name[256];
282     strcpy(name,argv[1]);
283     strcat(name,"_good");
284     DBRep::Set(name,aCompGood);
285   }
286   if (!aCompFailed.IsNull())
287   {
288     char name[256];
289     strcpy(name,argv[1]);
290     strcat(name,"_failed");
291     DBRep::Set(name,aCompFailed);
292   }
293   if (!aCompViolating.IsNull())
294   {
295     char name[256];
296     strcpy(name,argv[1]);
297     strcat(name,"_violating");
298     DBRep::Set(name,aCompViolating);
299   }
300
301   di<<"FastDiscret completed with :"<<"\n";
302   di<<"  MaxDeflection="<<maxdef<<"\n";
303   di<<"  NbNodes="<<nbnodes<<"\n";
304   di<<"  NbTriangles="<<nbtriangles<<"\n";
305   di<<"  NbFailed="<<nbfailed<<"\n";
306   di<<"  NbViolating="<<nbviolating<<"\n";
307
308   return 0;
309 }
310
311
312 //=======================================================================
313 //function : triangule
314 //purpose  : 
315 //=======================================================================
316
317
318 class BRepMesh_Couple
319 {
320 public:
321   BRepMesh_Couple() { myI1 = myI2 = 0; }
322   BRepMesh_Couple(const Standard_Integer I1,
323     const Standard_Integer I2)
324   { myI1 = I1; myI2 = I2; }
325
326   Standard_Integer myI1;
327   Standard_Integer myI2;
328 };
329
330 inline Standard_Boolean IsEqual(const BRepMesh_Couple& one,
331                                 const BRepMesh_Couple& other)
332 {
333   if (one.myI1 == other.myI1 &&
334     one.myI2 == other.myI2) return Standard_True;
335   else return Standard_False;
336 }
337
338 inline Standard_Integer HashCode(const BRepMesh_Couple& one,
339                                  const Standard_Integer Upper)
340 {
341   return ::HashCode((one.myI1+one.myI2), Upper);
342 }
343
344 typedef NCollection_Map<BRepMesh_Couple> BRepMesh_MapOfCouple;
345
346
347 static void AddLink(BRepMesh_MapOfCouple& aMap, 
348                     Standard_Integer v1,
349                     Standard_Integer v2)
350 {
351   Standard_Integer i1 = v1;
352   Standard_Integer i2 = v2;
353   if(i1 > i2) {
354     i1 = v2;
355     i2 = v1;
356   }
357   aMap.Add(BRepMesh_Couple(i1,i2));
358 }
359
360 static void MeshStats(const TopoDS_Shape& theSape,
361                       Standard_Integer& theNbTri,
362                       Standard_Integer& theNbEdges,
363                       Standard_Integer& theNbNodes)
364 {
365   theNbTri = 0;
366   theNbEdges = 0;
367   theNbNodes = 0;
368
369   Handle(Poly_Triangulation) T;
370   TopLoc_Location L;
371
372   for ( TopExp_Explorer ex(theSape, TopAbs_FACE); ex.More(); ex.Next()) {
373     TopoDS_Face F = TopoDS::Face(ex.Current());
374     T = BRep_Tool::Triangulation(F, L);
375     if (!T.IsNull()) {
376       theNbTri += T->NbTriangles();
377       theNbNodes += T->NbNodes();
378
379       BRepMesh_MapOfCouple aMap;
380       //count number of links
381       Poly_Array1OfTriangle& Trian = T->ChangeTriangles();
382       for(Standard_Integer i = 1; i<=Trian.Length();i++) {
383         Standard_Integer v1, v2, v3;
384         Trian(i).Get(v1,v2,v3);
385
386         AddLink(aMap, v1, v2);
387         AddLink(aMap, v2, v3);
388         AddLink(aMap, v3, v1);
389       }
390
391       theNbEdges+=aMap.Extent();
392     }
393   }
394 }
395
396 static Standard_Integer triangule(Draw_Interpretor& di, Standard_Integer nbarg, const char** argv)
397 {
398   if (nbarg < 4) return 1;
399
400   Standard_Boolean save = Standard_False;
401
402   const char *id1 = argv[2];
403   TopoDS_Shape S = DBRep::Get(id1);
404   if (S.IsNull()) return 1;
405   di << argv[1] << " ";
406   Standard_Real Deflect=Draw::Atof(argv[3]);
407   if (Deflect<=0.) {
408     di << " Donner la fleche !" << "\n";
409     return 1;
410   }
411
412   if (nbarg >4) {
413     save = (Draw::Atoi(argv[4])==1);
414   }
415
416   Standard_Boolean partage=Standard_True;
417   if (nbarg>5) {
418     partage=Draw::Atoi(argv[5])==1;
419   }
420
421   Handle(MeshTest_DrawableMesh) DM =
422     new MeshTest_DrawableMesh(S,Deflect,partage, save);
423
424   Draw::Set(argv[1],DM);
425
426   Standard_Integer nbn, nbl, nbe;
427   MeshStats(S, nbe, nbl, nbn);
428
429   di<<"(Resultat ("<<nbe<<" mailles) ("<<nbl<<" aretes) ("<<nbn<<" sommets))"<<"\n";
430
431   // passe de verification du maillage.
432   /*Standard_Integer nbc;
433   for (Standard_Integer iLi=1; iLi<= DM->Mesh()->NbEdges(); iLi++) {
434   const BRepMesh_Edge& ed=DM->Mesh()->Edge(iLi);
435   if (ed.Movability()!=BRepMesh_Deleted) {
436   nbc=struc->ElemConnectedTo(iLi).Extent();
437   if (nbc != 1 && nbc != 2) di <<"ERROR MAILLAGE Edge no "<< iLi<<"\n";
438   }
439   }*/
440
441   Bnd_Box bobo;
442
443   for (Standard_Integer lepnt=1; lepnt<DM->Mesh()->NbPoint3d(); lepnt++) {
444     bobo.Add(DM->Mesh()->Point3d(lepnt));
445   }
446   Standard_Real x,y,z,X,Y,Z;
447   bobo.Get(x,y,z,X,Y,Z);
448   Standard_Real delta=Max(X-x,Max(Y-y,Z-z));
449   if (delta>0) delta=Deflect/delta;
450   di << " Fleche de " << delta << " fois la taille de l''objet." << "\n";
451
452   return 0;
453 }
454
455 //=======================================================================
456 //function : addshape
457 //purpose  : 
458 //=======================================================================
459
460 Standard_Integer addshape(Draw_Interpretor&, Standard_Integer n, const char** a)
461 {
462   if (n < 3) return 1;
463   Handle(MeshTest_DrawableMesh) D =
464     Handle(MeshTest_DrawableMesh)::DownCast(Draw::Get(a[1]));
465   if (D.IsNull()) return 1;
466   TopoDS_Shape S = DBRep::Get(a[2]);
467   if (S.IsNull()) return 1;
468
469   D->Add(S);
470   Draw::Repaint();
471
472   return 0;
473 }
474
475
476 //=======================================================================
477 //function : smooth
478 //purpose  : 
479 //=======================================================================
480
481 /*Standard_Integer smooth(Draw_Interpretor&, Standard_Integer n, const char** a)
482 {
483 if (n < 2) return 1;
484 Handle(MeshTest_DrawableMesh) D =
485 Handle(MeshTest_DrawableMesh)::DownCast(Draw::Get(a[1]));
486 if (D.IsNull()) return 1;
487 Handle(BRepMesh_DataStructureOfDelaun) struc=
488 D->Mesh()->Result();
489 BRepMesh_Array1OfVertexOfDelaun toto(1,1);
490 BRepMesh_Delaun trial(struc, 
491 toto,
492 Standard_True);
493 trial.SmoothMesh(0.1);
494 Draw::Repaint();
495 return 0;
496 }
497 */
498
499 //=======================================================================
500 //function : edges
501 //purpose  : 
502 //=======================================================================
503
504 /*static Standard_Integer edges (Draw_Interpretor&, Standard_Integer n, const char** a)
505 {
506 if (n < 3) return 1;
507
508 Handle(MeshTest_DrawableMesh) D =
509 Handle(MeshTest_DrawableMesh)::DownCast(Draw::Get(a[1]));
510 if (D.IsNull()) return 1;
511 TopoDS_Shape S = DBRep::Get(a[2]);
512 if (S.IsNull()) return 1;
513
514 TopExp_Explorer ex;
515 TColStd_SequenceOfInteger& eseq = D->Edges();
516 Handle(BRepMesh_FastDiscret) M = D->Mesh();
517 Handle(BRepMesh_DataStructureOfDelaun) DS = M->Result();
518 Standard_Integer e1, e2, e3, iTri;
519 Standard_Boolean o1, o2, o3;
520
521 // the faces
522 for (ex.Init(S,TopAbs_FACE);ex.More();ex.Next()) {
523 const BRepMesh_MapOfInteger& elems = DS->ElemOfDomain();
524 BRepMesh_MapOfInteger::Iterator it;
525 for (it.Initialize(elems); it.More(); it.Next()) {
526 iTri = it.Key();
527 const BRepMesh_Triangle& triang = M->Triangle(iTri);
528 if (triang.Movability()!=BRepMesh_Deleted) {
529 triang.Edges(e1, e2, e3, o1, o2, o3);
530 eseq.Append(e1);
531 eseq.Append(e2);
532 eseq.Append(e3);
533 }
534 }
535 }
536
537 // the edges
538 //for (ex.Init(S,TopAbs_EDGE,TopAbs_FACE);ex.More();ex.Next()) {
539 //}
540
541 Draw::Repaint();
542 return 0;
543 }
544 */
545
546 //=======================================================================
547 //function : vertices
548 //purpose  : 
549 //=======================================================================
550
551 static Standard_Integer vertices (Draw_Interpretor&, Standard_Integer n, const char** a)
552 {
553   if (n < 3) return 1;
554
555   Handle(MeshTest_DrawableMesh) D =
556     Handle(MeshTest_DrawableMesh)::DownCast(Draw::Get(a[1]));
557   if (D.IsNull()) return 1;
558   TopoDS_Shape S = DBRep::Get(a[2]);
559   if (S.IsNull()) return 1;
560
561   TopExp_Explorer ex;
562   TColStd_SequenceOfInteger& vseq = D->Vertices();
563   Handle(BRepMesh_FastDiscret) M = D->Mesh();
564
565   // the faces
566   for (ex.Init(S,TopAbs_FACE);ex.More();ex.Next()) {
567     BRepMesh_MapOfInteger vtx;
568     M->VerticesOfDomain(vtx);
569     for (BRepMesh_MapOfInteger::Iterator it(vtx); it.More(); it.Next())
570       vseq.Append(it.Key());
571   }
572
573
574   // the edges
575   //for (ex.Init(S,TopAbs_EDGE,TopAbs_FACE);ex.More();ex.Next()) {
576   //}
577
578   Draw::Repaint();
579   return 0;
580 }
581
582 //=======================================================================
583 //function : medge
584 //purpose  : 
585 //=======================================================================
586
587 static Standard_Integer medge (Draw_Interpretor&, Standard_Integer n, const char** a)
588 {
589   if (n < 3) return 1;
590
591   Handle(MeshTest_DrawableMesh) D =
592     Handle(MeshTest_DrawableMesh)::DownCast(Draw::Get(a[1]));
593   if (D.IsNull()) return 1;
594
595   Standard_Integer i,j,e;
596   TColStd_SequenceOfInteger& eseq = D->Edges();
597   for (i = 2; i < n; i++) {
598     e = Draw::Atoi(a[i]);
599     if (e > 0)
600       eseq.Append(e);
601     else if (e < 0) {
602       e = -e;
603       j = 1; 
604       while (j <= eseq.Length()) {
605         if (eseq(j) == e) 
606           eseq.Remove(j);
607         else
608           j++;
609       }
610     }
611     else
612       eseq.Clear();
613   }
614
615   Draw::Repaint();
616   return 0;
617 }
618
619
620 //=======================================================================
621 //function : mvertex
622 //purpose  : 
623 //=======================================================================
624
625 static Standard_Integer mvertex (Draw_Interpretor&, Standard_Integer n, const char** a)
626 {
627   if (n < 3) return 1;
628
629   Handle(MeshTest_DrawableMesh) D =
630     Handle(MeshTest_DrawableMesh)::DownCast(Draw::Get(a[1]));
631   if (D.IsNull()) return 1;
632
633   Standard_Integer i,j,v;
634   TColStd_SequenceOfInteger& vseq = D->Vertices();
635   for (i = 2; i < n; i++) {
636     v = Draw::Atoi(a[i]);
637     if (v > 0)
638       vseq.Append(v);
639     else if (v < 0) {
640       v = -v;
641       j = 1;
642       while (j <= vseq.Length()) {
643         if (vseq(j) == v)
644           vseq.Remove(v);
645         else
646           j++;
647       }
648     }
649     else
650       vseq.Clear();
651   }
652   Draw::Repaint();
653   return 0;
654 }
655
656
657 //=======================================================================
658 //function : triangle
659 //purpose  : 
660 //=======================================================================
661
662 static Standard_Integer triangle (Draw_Interpretor&, Standard_Integer n, const char** a)
663 {
664   if (n < 3) return 1;
665
666   Handle(MeshTest_DrawableMesh) D =
667     Handle(MeshTest_DrawableMesh)::DownCast(Draw::Get(a[1]));
668   if (D.IsNull()) return 1;
669
670   Standard_Integer i,j,v;
671   TColStd_SequenceOfInteger& tseq = D->Triangles();
672   for (i = 2; i < n; i++) {
673     v = Draw::Atoi(a[i]);
674     if (v > 0)
675       tseq.Append(v);
676     else if (v < 0) {
677       v = -v;
678       j = 1;
679       while (j <= tseq.Length()) {
680         if (tseq(j) == v)
681           tseq.Remove(v);
682         else
683           j++;
684       }
685     }
686     else
687       tseq.Clear();
688   }
689   Draw::Repaint();
690   return 0;
691 }
692
693 //=======================================================================
694 //function : dumpvertex
695 //purpose  : 
696 //=======================================================================
697
698 /*
699 Standard_Integer dumpvertex(Draw_Interpretor& di, Standard_Integer argc, const char** argv)
700 {
701 if (argc < 2) return 1;
702
703 Handle(MeshTest_DrawableMesh) D =
704 Handle(MeshTest_DrawableMesh)::DownCast(Draw::Get(argv[1]));
705 if (D.IsNull()) return 1;
706
707 Handle(BRepMesh_DataStructureOfDelaun) struc = D->Mesh()->Result();
708
709 Standard_Integer in=1;
710 if (argc>=3) {
711 in=Draw::Atoi(argv[2]);
712 in=Max(1,in);
713 }
714 Standard_Integer nbn=in;
715 if (argc>=4) {
716 nbn=Draw::Atoi(argv[3]);
717 nbn=Min(nbn,struc->NbNodes());
718 }
719
720 for (; in<=nbn; in++) {
721 BRepMesh_Vertex nod=struc->GetNode(in);
722 di<<"(node "<<in<<" (uv "<<nod.Coord().X()
723 <<" "<<nod.Coord().Y()<<") (3d "
724 <<nod.Location3d()<<") ";
725 printdegree(nod.Movability(), di);
726 di<<" (edgeconex";
727 BRepMesh_ListOfInteger::Iterator tati(struc->LinkNeighboursOf(in));
728 for (; tati.More(); tati.Next()) di<<" "<<tati.Value();
729 di << "))\n";
730 }
731 di <<"\n";
732 return 0;
733 }
734
735 //=======================================================================
736 //function : dumpedge
737 //purpose  : 
738 //=======================================================================
739
740 Standard_Integer dumpedge(Draw_Interpretor& di, Standard_Integer argc, const char** argv)
741 {
742 if (argc < 2) return 1;
743
744 Handle(MeshTest_DrawableMesh) D =
745 Handle(MeshTest_DrawableMesh)::DownCast(Draw::Get(argv[1]));
746 if (D.IsNull()) return 1;
747
748 Handle(BRepMesh_DataStructureOfDelaun) struc=D->Mesh()->Result();
749 Standard_Integer il=1;
750 if (argc>=3) {
751 il=Draw::Atoi(argv[2]);
752 il=Max(1, il);
753 }
754 Standard_Integer nbl=il;
755 if (argc>=4) {
756 nbl=Draw::Atoi(argv[3]);
757 nbl=Min(nbl, struc->NbLinks());
758 }
759
760 for (; il<=nbl; il++) {
761 BRepMesh_Edge edg=struc->GetLink(il);
762 di << "(edge "<<il<<" ("<<edg.FirstNode()<<" "<<edg.LastNode()
763 <<" ";
764 printdegree(edg.Movability(), di);
765 di<<") (triconex";
766 const BRepMesh_PairOfIndex& pair = struc->ElemConnectedTo(il);
767 for (Standard_Integer j = 1, jn = pair.Extent(); j <= jn; j++)
768 di<<" "<<pair.Index(j);
769 di << "))\n";
770 }
771 di <<"\n";
772 return 0;
773 }
774
775 //=======================================================================
776 //function : dumptriangle
777 //purpose  : 
778 //=======================================================================
779
780 Standard_Integer dumptriangle(Draw_Interpretor& di, Standard_Integer argc, const char** argv)
781 {
782 if (argc < 2) return 1;
783
784 Handle(MeshTest_DrawableMesh) D =
785 Handle(MeshTest_DrawableMesh)::DownCast(Draw::Get(argv[1]));
786 if (D.IsNull()) return 1;
787
788 Handle(BRepMesh_DataStructureOfDelaun) struc=D->Mesh()->Result();
789 Standard_Integer ie=1;
790 if (argc>=3) {
791 ie=Draw::Atoi(argv[2]);
792 ie=Max(1, ie);
793 }
794 Standard_Integer nbe=ie;
795 if (argc>=4) {
796 nbe=Draw::Atoi(argv[3]);
797 nbe=Min(nbe, struc->NbElements());
798 }
799
800 Standard_Integer e1, e2, e3;
801 Standard_Boolean o1, o2, o3;
802
803 for (; ie<=nbe; ie++) {
804 BRepMesh_Triangle tri=struc->GetElement(ie);
805 tri.Edges(e1, e2, e3, o1, o2, o3); 
806 if (o1) e1=-e1;
807 if (o2) e2=-e2;
808 if (o3) e3=-e3;
809 di<<" (maille "<<ie<<" (links "<<e1<<" "
810 <<e2<<" "<<e3<<")";
811 printdegree(tri.Movability(), di);
812 di<<")\n";
813 }
814 di << "\n";
815 return 0;
816 }
817 */
818
819 //=======================================================================
820 //function : trianglesinfo
821 //purpose  : 
822 //=======================================================================
823 static Standard_Integer trianglesinfo(Draw_Interpretor& di, Standard_Integer n, const char** a)
824 {
825   if (n != 2) return 1;
826   TopoDS_Shape S = DBRep::Get(a[1]);
827   if (S.IsNull()) return 1;
828   TopExp_Explorer ex;
829   Handle(Poly_Triangulation) T;
830   TopLoc_Location L;
831
832   Standard_Real MaxDeflection = 0.0;
833   Standard_Integer nbtriangles = 0, nbnodes = 0;
834   for (ex.Init(S, TopAbs_FACE); ex.More(); ex.Next()) {
835     TopoDS_Face F = TopoDS::Face(ex.Current());
836     T = BRep_Tool::Triangulation(F, L);
837     if (!T.IsNull()) {
838       nbtriangles += T->NbTriangles();
839       nbnodes += T->NbNodes();
840       if (T->Deflection() > MaxDeflection)
841         MaxDeflection = T->Deflection();
842     }
843   }
844
845   di<<"\n";
846   di<<"This shape contains " <<nbtriangles<<" triangles."<<"\n";
847   di<<"                    " <<nbnodes    <<" nodes."<<"\n";
848   di<<"Maximal deflection " <<MaxDeflection<<"\n";
849   di<<"\n";
850 #ifdef DEB_MESH_CHRONO
851   Standard_Real tot, addp, unif, contr, inter;
852   Standard_Real edges, mailledges, etuinter, lastcontrol, stock;
853   Standard_Real add11, add12, add2, upda, pointvalid;
854   Standard_Real isos, pointsisos;
855   chTotal.Show(tot); chAddPoint.Show(addp); chUnif.Show(unif); 
856   chControl.Show(contr); chInternal.Show(inter);
857   chEdges.Show(edges); chMaillEdges.Show(mailledges);
858   chEtuInter.Show(etuinter); chLastControl.Show(lastcontrol); 
859   chStock.Show(stock);
860   chAdd11.Show(add11); chAdd12.Show(add12); chAdd2.Show(add2); chUpdate.Show(upda);
861   chPointValid.Show(pointvalid); chIsos.Show(isos); chPointsOnIsos.Show(pointsisos);
862
863   if (tot > 0.00001) {
864     di <<"temps total de maillage:     "<<tot        <<" seconds"<< "\n";
865     di <<"dont: "<< "\n";
866     di <<"discretisation des edges:    "<<edges      <<" seconds---> "<< 100*edges/tot      <<" %"<<"\n";
867     di <<"maillage des edges:          "<<mailledges <<" seconds---> "<< 100*mailledges/tot <<" %"<<"\n";
868     di <<"controle et points internes: "<<etuinter   <<" seconds---> "<< 100*etuinter/tot   <<" %"<<"\n";
869     di <<"derniers controles:          "<<lastcontrol<<" seconds---> "<< 100*lastcontrol/tot<<" %"<<"\n";
870     di <<"stockage dans la S.D.        "<<stock      <<" seconds---> "<< 100*stock/tot      <<" %"<<"\n";
871     di << "\n";
872     di <<"et plus precisement: "<<"\n";
873     di <<"Add 11ere partie :           "<<add11     <<" seconds---> "<<100*add11/tot      <<" %"<<"\n";
874     di <<"Add 12ere partie :           "<<add12     <<" seconds---> "<<100*add12/tot      <<" %"<<"\n";
875     di <<"Add 2eme partie :            "<<add2      <<" seconds---> "<<100*add2/tot       <<" %"<<"\n";
876     di <<"Update :                     "<<upda      <<" seconds---> "<<100*upda/tot       <<" %"<<"\n";
877     di <<"AddPoint :                   "<<addp      <<" seconds---> "<<100*addp/tot       <<" %"<<"\n";
878     di <<"UniformDeflection            "<<unif      <<" seconds---> "<<100*unif/tot       <<" %"<<"\n";
879     di <<"Controle :                   "<<contr     <<" seconds---> "<<100*contr/tot      <<" %"<<"\n";
880     di <<"Points Internes:             "<<inter     <<" seconds---> "<<100*inter/tot      <<" %"<<"\n";
881     di <<"calcul des isos et du, dv:   "<<isos      <<" seconds---> "<<100*isos/tot       <<" %"<<"\n";
882     di <<"calcul des points sur isos:  "<<pointsisos<<" seconds---> "<<100*pointsisos/tot <<" %"<<"\n";
883     di <<"IsPointValid:                "<<pointvalid<<" seconds---> "<<100*pointvalid/tot <<" %"<<"\n";
884     di << "\n";
885
886
887     di <<"nombre d'appels de controle apres points internes          : "<< NbControls << "\n";
888     di <<"nombre de points sur restrictions                          : "<< D0Edges    << "\n";
889     di <<"nombre de points calcules par UniformDeflection            : "<< D0Unif     << "\n";
890     di <<"nombre de points calcules dans InternalVertices            : "<< D0Internal << "\n";
891     di <<"nombre de points calcules dans Control                     : "<< D0Control  << "\n";
892     if (nbnodes-D0Edges != 0) { 
893       Standard_Real ratio = (Standard_Real)(D0Internal+D0Control)/ (Standard_Real)(nbnodes-D0Edges);
894       di <<"---> Ratio: (D0Internal+D0Control) / (nbNodes-nbOnEdges)   : "<< ratio      << "\n";
895     }
896
897     di << "\n";
898
899     chTotal.Reset(); chAddPoint.Reset(); chUnif.Reset(); 
900     chControl.Reset(); chInternal.Reset();
901     chEdges.Reset(); chMaillEdges.Reset();
902     chEtuInter.Reset(); chLastControl.Reset(); 
903     chStock.Reset();
904     chAdd11.Reset(); chAdd12.Reset(); chAdd2.Reset(); chUpdate.Reset();
905     chPointValid.Reset(); chIsos.Reset(); chPointsOnIsos.Reset();
906
907   }
908 #endif
909   return 0;
910 }
911
912 //=======================================================================
913 //function : veriftriangles
914 //purpose  : 
915 //=======================================================================
916
917 static Standard_Integer veriftriangles(Draw_Interpretor& di, Standard_Integer n, const char** a)
918 {
919   if (n < 2) return 1;
920   Standard_Boolean quiet = 1;
921   if (n == 3) quiet = 0;
922   TopoDS_Shape Sh = DBRep::Get(a[1]);
923   if (Sh.IsNull()) return 1;
924   TopExp_Explorer ex;
925   Handle(Poly_Triangulation) T;
926   TopLoc_Location L;
927   Standard_Integer i, n1, n2, n3;
928   gp_Pnt2d mitri, v1, v2, v3, mi2d1, mi2d2, mi2d3;
929   gp_XYZ vecEd1, vecEd2, vecEd3;
930   //  Standard_Real dipo, dm, dv, d1, d2, d3, defle;
931   Standard_Real dipo, dv, d1, d2, d3, defle;
932   Handle(Geom_Surface) S;
933   Standard_Integer nbface = 0;
934   gp_Pnt PP;
935
936   for (ex.Init(Sh, TopAbs_FACE); ex.More(); ex.Next()) {
937     TopoDS_Face F = TopoDS::Face(ex.Current());
938     nbface++;
939     T = BRep_Tool::Triangulation(F, L);
940     Standard_Real deflemax = 0, deflemin = 1.e100;
941     if (!T.IsNull()) {
942       Standard_Real defstock = T->Deflection();
943       const Poly_Array1OfTriangle& triangles  = T->Triangles();
944       const TColgp_Array1OfPnt2d&  Nodes2d    = T->UVNodes();
945       const TColgp_Array1OfPnt&    Nodes      = T->Nodes();
946
947       S = BRep_Tool::Surface(F, L);
948
949       for(i = 1; i <= triangles.Length(); i++) {
950         if (F.Orientation() == TopAbs_REVERSED) 
951           triangles(i).Get(n1,n3,n2);
952         else 
953           triangles(i).Get(n1,n2,n3);
954
955         const gp_XY& xy1 = Nodes2d(n1).XY();
956         const gp_XY& xy2 = Nodes2d(n2).XY();
957         const gp_XY& xy3 = Nodes2d(n3).XY();
958
959         mi2d1.SetCoord((xy2.X()+xy3.X())*0.5, 
960           (xy2.Y()+xy3.Y())*0.5);
961         mi2d2.SetCoord((xy1.X()+xy3.X())*0.5, 
962           (xy1.Y()+xy3.Y())*0.5);
963         mi2d3.SetCoord((xy1.X()+xy2.X())*0.5, 
964           (xy1.Y()+xy2.Y())*0.5);
965
966         gp_XYZ p1 = Nodes(n1).Transformed(L.Transformation()).XYZ();
967         gp_XYZ p2 = Nodes(n2).Transformed(L.Transformation()).XYZ();
968         gp_XYZ p3 = Nodes(n3).Transformed(L.Transformation()).XYZ();
969
970         vecEd1=p2-p1;
971         vecEd2=p3-p2;
972         vecEd3=p1-p3;
973         d1=vecEd1.SquareModulus();
974         d2=vecEd2.SquareModulus();
975         d3=vecEd3.SquareModulus();
976
977         if (d1!=0. && d2!=0. && d3!=0.) {
978           gp_XYZ equa(vecEd1^vecEd2);
979           dv=equa.Modulus();
980           if (dv>0.) {
981             equa.SetCoord(equa.X()/dv, equa.Y()/dv, equa.Z()/dv);
982             dipo=equa*p1;
983
984
985             mitri.SetCoord(ONETHIRD*(xy1.X()+xy2.X()+xy3.X()),
986               ONETHIRD*(xy1.Y()+xy2.Y()+xy3.Y()));
987             v1.SetCoord(ONETHIRD*mi2d1.X()+TWOTHIRD*xy1.X(), 
988               ONETHIRD*mi2d1.Y()+TWOTHIRD*xy1.Y());
989             v2.SetCoord(ONETHIRD*mi2d2.X()+TWOTHIRD*xy2.X(), 
990               ONETHIRD*mi2d2.Y()+TWOTHIRD*xy2.Y());
991             v3.SetCoord(ONETHIRD*mi2d3.X()+TWOTHIRD*xy3.X(), 
992               ONETHIRD*mi2d3.Y()+TWOTHIRD*xy3.Y());
993
994             S->D0(mi2d1.X(), mi2d1.Y(), PP);
995             PP = PP.Transformed(L.Transformation());
996             defle = Abs((equa*PP.XYZ())-dipo);
997             deflemax = Max(deflemax, defle);
998             deflemin = Min(deflemin, defle);
999
1000             S->D0(mi2d2.X(), mi2d2.Y(), PP);
1001             PP = PP.Transformed(L.Transformation());
1002             defle = Abs((equa*PP.XYZ())-dipo);
1003             deflemax = Max(deflemax, defle);
1004             deflemin = Min(deflemin, defle);
1005
1006             S->D0(mi2d3.X(), mi2d3.Y(), PP);
1007             PP = PP.Transformed(L.Transformation());
1008             defle = Abs((equa*PP.XYZ())-dipo);
1009             deflemax = Max(deflemax, defle);
1010             deflemin = Min(deflemin, defle);
1011
1012             S->D0(v1.X(), v1.Y(), PP);
1013             PP = PP.Transformed(L.Transformation());
1014             defle = Abs((equa*PP.XYZ())-dipo);
1015             deflemax = Max(deflemax, defle);
1016             deflemin = Min(deflemin, defle);
1017
1018             S->D0(v2.X(), v2.Y(), PP);
1019             PP = PP.Transformed(L.Transformation());
1020             defle = Abs((equa*PP.XYZ())-dipo);
1021             deflemax = Max(deflemax, defle);
1022             deflemin = Min(deflemin, defle);
1023
1024             S->D0(v3.X(), v3.Y(), PP);
1025             PP = PP.Transformed(L.Transformation());
1026             defle = Abs((equa*PP.XYZ())-dipo);
1027             deflemax = Max(deflemax, defle);
1028             deflemin = Min(deflemin, defle);
1029
1030             S->D0(mitri.X(), mitri.Y(), PP);
1031             PP = PP.Transformed(L.Transformation());
1032             defle = Abs((equa*PP.XYZ())-dipo);
1033             deflemax = Max(deflemax, defle);
1034             deflemin = Min(deflemin, defle);
1035
1036             if (defle > defstock) {
1037               di <<"face "<< nbface <<" deflection = " << defle <<" pour "<<defstock <<" stockee."<<"\n";
1038             }
1039           }
1040         }
1041       }
1042       if (!quiet) {
1043         di <<"face "<< nbface<<", deflemin = "<< deflemin<<", deflemax = "<<deflemax<<"\n";
1044       }
1045
1046     }
1047   }
1048
1049
1050   return 0;
1051 }
1052
1053
1054
1055
1056 //=======================================================================
1057 //function : tri2d
1058 //purpose  : 
1059 //=======================================================================
1060
1061 Standard_Integer tri2d(Draw_Interpretor&, Standard_Integer n, const char** a)
1062 {
1063
1064   if (n != 2) return 1;
1065   TopoDS_Shape aLocalShape = DBRep::Get(a[1]);
1066   TopoDS_Face F = TopoDS::Face(aLocalShape);
1067   //  TopoDS_Face F = TopoDS::Face(DBRep::Get(a[1]));
1068   if (F.IsNull()) return 1;
1069   Handle(Poly_Triangulation) T;
1070   TopLoc_Location L;
1071
1072   T = BRep_Tool::Triangulation(F, L);
1073   if (!T.IsNull()) {
1074     // Build the connect tool
1075     Poly_Connect pc(T);
1076
1077     Standard_Integer i,j, nFree, nInternal, nbTriangles = T->NbTriangles();
1078     Standard_Integer t[3];
1079
1080     // count the free edges
1081     nFree = 0;
1082     for (i = 1; i <= nbTriangles; i++) {
1083       pc.Triangles(i,t[0],t[1],t[2]);
1084       for (j = 0; j < 3; j++)
1085         if (t[j] == 0) nFree++;
1086     }
1087
1088     // allocate the arrays
1089     TColStd_Array1OfInteger Free(1,2*nFree);
1090     nInternal = (3*nbTriangles - nFree) / 2;
1091     TColStd_Array1OfInteger Internal(0,2*nInternal);
1092
1093     Standard_Integer fr = 1, in = 1;
1094     const Poly_Array1OfTriangle& triangles = T->Triangles();
1095     Standard_Integer nodes[3];
1096     for (i = 1; i <= nbTriangles; i++) {
1097       pc.Triangles(i,t[0],t[1],t[2]);
1098       triangles(i).Get(nodes[0],nodes[1],nodes[2]);
1099       for (j = 0; j < 3; j++) {
1100         Standard_Integer k = (j+1) % 3;
1101         if (t[j] == 0) {
1102           Free(fr)   = nodes[j];
1103           Free(fr+1) = nodes[k];
1104           fr += 2;
1105         }
1106         // internal edge if this triangle has a lower index than the adjacent
1107         else if (i < t[j]) {
1108           Internal(in)   = nodes[j];
1109           Internal(in+1) = nodes[k];
1110           in += 2;
1111         }
1112       }
1113     }
1114
1115     // Display the edges
1116     if (T->HasUVNodes()) {
1117       const TColgp_Array1OfPnt2d& Nodes2d = T->UVNodes();
1118
1119       Handle(Draw_Segment2D) Seg;
1120
1121       // free edges
1122       Standard_Integer nn;
1123       nn = Free.Length() / 2;
1124       for (i = 1; i <= nn; i++) {
1125         Seg = new Draw_Segment2D(Nodes2d(Free(2*i-1)),
1126           Nodes2d(Free(2*i)),
1127           Draw_rouge);
1128         dout << Seg;
1129       }
1130
1131       // internal edges
1132
1133       nn = nInternal;
1134       for (i = 1; i <= nn; i++) {
1135         Seg = new Draw_Segment2D(Nodes2d(Internal(2*i-1)),
1136           Nodes2d(Internal(2*i)),
1137           Draw_bleu);
1138         dout << Seg;
1139       }
1140     }
1141     dout.Flush();
1142   }
1143
1144   return 0;
1145 }
1146
1147
1148
1149
1150 //=======================================================================
1151 //function : wavefront
1152 //purpose  : 
1153 //=======================================================================
1154
1155 static Standard_Integer wavefront(Draw_Interpretor&, Standard_Integer nbarg, const char** argv)
1156 {
1157   if (nbarg < 2) return 1;
1158
1159   TopoDS_Shape S = DBRep::Get(argv[1]);
1160   if (S.IsNull()) return 1;
1161
1162   // creation du maillage s'il n'existe pas.
1163
1164   Bnd_Box B;
1165   Standard_Real aXmin, aYmin, aZmin, aXmax, aYmax, aZmax;
1166   BRepBndLib::Add(S, B);
1167   B.Get(aXmin, aYmin, aZmin, aXmax, aYmax, aZmax);
1168   Standard_Real aDeflection = 
1169     MAX3( aXmax-aXmin , aYmax-aYmin , aZmax-aZmin) * 0.004;
1170
1171   BRepMesh_IncrementalMesh aMesh (S, aDeflection);
1172
1173
1174   TopLoc_Location L;
1175   TopExp_Explorer ex;
1176
1177   Standard_Integer i, nbface = 0;
1178   Standard_Boolean OK = Standard_True;
1179   gp_Vec D1U,D1V;
1180   gp_Vec D2U,D2V,D2UV;
1181   gp_Dir Nor;
1182   gp_Pnt P;
1183   Standard_Real U, V;
1184   CSLib_DerivativeStatus Status;
1185   CSLib_NormalStatus NStat;
1186   Standard_Real x, y, z;
1187   Standard_Integer n1, n2, n3;
1188   Standard_Integer k1, k2, k3;
1189
1190   char ffile[100];
1191
1192   if (nbarg == 3) {
1193     strcpy(ffile, argv[2]);
1194     strcat(ffile, ".obj");
1195   }
1196   else  strcpy(ffile, "wave.obj");
1197   FILE* outfile = fopen(ffile, "w");
1198
1199
1200   fprintf(outfile, "%s  %s\n%s %s\n\n", "# CASCADE   ","MATRA DATAVISION", "#", ffile);
1201
1202   Standard_Integer nbNodes, totalnodes = 0, nbpolygons = 0;
1203   for (ex.Init(S, TopAbs_FACE); ex.More(); ex.Next()) {
1204     nbface++;
1205     TopoDS_Face F = TopoDS::Face(ex.Current());
1206     Handle(Poly_Triangulation) Tr = BRep_Tool::Triangulation(F, L);
1207
1208     if (!Tr.IsNull()) {
1209       nbNodes = Tr->NbNodes();
1210       const TColgp_Array1OfPnt& Nodes = Tr->Nodes();
1211
1212       // les noeuds.
1213       for (i = 1; i <= nbNodes; i++) {
1214         gp_Pnt Pnt = Nodes(i).Transformed(L.Transformation());
1215         x = Pnt.X();
1216         y = Pnt.Y();
1217         z = Pnt.Z();
1218         fprintf(outfile, "%s      %f  %f  %f\n", "v", x, y, z);
1219       }
1220
1221       fprintf(outfile, "\n%s    %d\n\n", "# number of vertex", nbNodes);
1222
1223
1224       // les normales.
1225
1226       if (Tr->HasUVNodes()) {
1227         const TColgp_Array1OfPnt2d& UVNodes = Tr->UVNodes();
1228         BRepAdaptor_Surface BS(F, Standard_False);
1229
1230         for (i = 1; i <= nbNodes; i++) {
1231           U = UVNodes(i).X();
1232           V = UVNodes(i).Y();
1233
1234           BS.D1(U,V,P,D1U,D1V);
1235           CSLib::Normal(D1U,D1V,Precision::Angular(),Status,Nor);
1236           if (Status != CSLib_Done) {
1237             BS.D2(U,V,P,D1U,D1V,D2U,D2V,D2UV);
1238             CSLib::Normal(D1U,D1V,D2U,D2V,D2UV,Precision::Angular(),OK,NStat,Nor);
1239           }
1240           if (F.Orientation() == TopAbs_REVERSED) Nor.Reverse();
1241
1242           fprintf(outfile, "%s      %f  %f  %f\n", "vn", Nor.X(), Nor.Y(), Nor.Z());
1243         }
1244
1245         fprintf(outfile, "\n%s    %d\n\n", "# number of vertex normals", nbNodes);
1246       }
1247
1248       fprintf(outfile, "%s    %d\n", "s", nbface);
1249
1250       // les triangles.
1251       Standard_Integer nbTriangles = Tr->NbTriangles();
1252       const Poly_Array1OfTriangle& triangles = Tr->Triangles();
1253
1254
1255       for (i = 1; i <= nbTriangles; i++) {
1256         if (F.Orientation()  == TopAbs_REVERSED)
1257           triangles(i).Get(n1, n3, n2);
1258         else 
1259           triangles(i).Get(n1, n2, n3);
1260         k1 = n1+totalnodes;
1261         k2 = n2+totalnodes;
1262         k3 = n3+totalnodes;
1263         fprintf(outfile, "%s %d%s%d %d%s%d %d%s%d\n", "fo", k1,"//", k1, k2,"//", k2, k3,"//", k3);
1264       }
1265       nbpolygons += nbTriangles;
1266       totalnodes += nbNodes;
1267
1268       fprintf(outfile, "\n%s    %d\n", "# number of smooth groups", nbface);
1269       fprintf(outfile, "\n%s    %d\n", "# number of polygons", nbpolygons);
1270
1271     }
1272   }
1273
1274   fclose(outfile);
1275
1276   return 0;
1277 }
1278
1279
1280 //=======================================================================
1281 //function : onetriangulation
1282 //purpose  : 
1283 //=======================================================================
1284
1285 Standard_Integer onetriangulation(Draw_Interpretor&, Standard_Integer /*nbarg*/, const char** /*argv*/)
1286 {
1287
1288   /*
1289
1290   if (nbarg < 2) return 1;
1291
1292   TopoDS_Shape S = DBRep::Get(argv[1]);
1293   if (S.IsNull()) return 1;
1294
1295   Handle(Poly_Triangulation) TFinale;
1296   char name[100];
1297   Standard_Integer nbshell = 0;
1298
1299   TopExp_Explorer ex, exs, ex2;
1300
1301   for (ex.Init(S, TopAbs_SHELL); ex.More(); ex.Next()) {
1302   nbshell++;
1303   TopoDS_Shell Sh = TopoDS::Shell(ex.Current());
1304
1305   for (exs.Init(Sh, TopAbs_Face); exs.More(); exs.Next()) {
1306   TopoDS_Face F = TopoDS::Face(exs.Current());
1307   Handle(Poly_Triangulation) T = BRep_Tool::Triangulation(F, L);
1308
1309   for (ex2.Init(F, TopAbs_EDGE); ex2.More(); ex2.Next()) {
1310   TopoDS_Edge edge = TopoDS::Edge(ex2.Current());
1311   const TColgp_Array1OfPnt& Nodes = T->Nodes();
1312   const Poly_Array1OfTriangle& triangles = T->Triangles();
1313
1314   if (mapedges.IsBound(edge)) {
1315   const TColStd_ListOfTransient& L = edges.Find(edge);
1316   const Handle(Poly_PolygonOnTriangulation)& P = 
1317   *(Handle(Poly_PolygonOnTriangulation)*)&(L.First());
1318   const TColStd_Array1OfInteger& NOD = P->Nodes();
1319
1320   }
1321   }
1322   }
1323
1324   Sprintf(name, "%s_%i", "tr", nbshell);
1325   DrawTrSurf::Set(name, TFinale);
1326
1327   }
1328
1329   */
1330   return 0;
1331 }
1332
1333
1334 #if 0
1335
1336 //=======================================================================
1337 //function : vb
1338 //purpose  : 
1339 //=======================================================================
1340
1341 Standard_Integer vb(Draw_Interpretor& di, Standard_Integer nbarg, const char** argv)
1342 {
1343   Standard_Integer NbPoints = 1, Deg = 1;
1344
1345   for (Deg = 1; Deg <= 25; Deg++) {
1346     for (NbPoints = 1; NbPoints <= 24; NbPoints++) {
1347
1348       math_Vector GaussP(1, NbPoints), GaussW(1, NbPoints);
1349       math_Vector TheWeights(1, NbPoints), VBParam(1, NbPoints);
1350       math_Matrix VB(1, Deg+1, 1, NbPoints);
1351
1352       math::GaussPoints(NbPoints, GaussP);
1353
1354       Standard_Integer i, j, classe = Deg+1, cl1 = Deg;
1355
1356       // calcul et mise en ordre des parametres et des poids:
1357       for (i = 1; i <= NbPoints; i++) {
1358         if (i <=  (NbPoints+1)/2) {
1359           VBParam(NbPoints-i+1)  = 0.5*(1 + GaussP(i));
1360         }
1361         else {
1362           VBParam(i-(NbPoints+1)/2)  = 0.5*(1 + GaussP(i));
1363         }
1364       }
1365
1366
1367       // Calcul du VB (Valeur des fonctions de Bernstein):
1368       for (i = 1; i <= classe; i++) {
1369         for (j = 1; j <= NbPoints; j++) {
1370           VB(i,j)=PLib::Binomial(cl1,i-1)*Pow((1-VBParam(j)),classe-i)*Pow(VBParam(j),i-1);
1371         }
1372       }
1373
1374
1375       for (i = 1; i <= classe; i++) {
1376         for (j = 1; j <= NbPoints; j++) {
1377           di<< VB(i, j) << ", ";
1378         }
1379       }
1380       di << "\n" << "\n";
1381     }
1382   }
1383   return 0;
1384 }  
1385 //=======================================================================
1386 //function : extrema
1387 //purpose  : 
1388 //=======================================================================
1389
1390 Standard_Integer extrema(Draw_Interpretor& di, Standard_Integer nbarg, const char** argv)
1391 {
1392
1393
1394   Handle(Geom_Curve) C = DrawTrSurf::GetCurve(argv[1]);
1395
1396   Standard_Real X, Y, Z, U0;
1397   X = Draw::Atof(argv[2]);
1398   Y = Draw::Atof(argv[3]);
1399   Z = Draw::Atof(argv[4]);
1400   U0 = Draw::Atof(argv[5]);
1401
1402   gp_Pnt P(X, Y, Z);
1403   GeomAdaptor_Curve GC(C);
1404   Standard_Real tol = 1.e-09;
1405   Extrema_LocateExtPC ext(P, GC, U0, tol);
1406
1407   if (ext.IsDone()) {
1408     gp_Pnt P1 = ext.Point().Value();
1409     di <<"distance =  "<<ext.Value() << "\n";
1410     di <<"point =     "<<P1.X()<<" "<<P1.Y()<<" "<< P1.Z()<< "\n";
1411     di <<"parametre = "<<ext.Point().Parameter()<<"\n";
1412   }
1413
1414   return 0;
1415 }
1416
1417 #endif
1418
1419
1420 //=======================================================================
1421 //function : triedgepoints
1422 //purpose  : 
1423 //=======================================================================
1424
1425 Standard_Integer triedgepoints(Draw_Interpretor& di, Standard_Integer nbarg, const char** argv)
1426 {
1427   if( nbarg < 2 )
1428     return 1;
1429
1430   for( Standard_Integer i = 1; i < nbarg; i++ )
1431   {
1432     TopoDS_Shape aShape = DBRep::Get(argv[i]);
1433     if ( aShape.IsNull() )
1434       continue;
1435
1436     Handle(Poly_PolygonOnTriangulation) aPoly;
1437     Handle(Poly_Triangulation)          aT;
1438     TopLoc_Location                     aLoc;
1439     TopTools_MapOfShape                 anEdgeMap;
1440     TopTools_MapIteratorOfMapOfShape    it;
1441     
1442     if( aShape.ShapeType() == TopAbs_EDGE )
1443     {
1444       anEdgeMap.Add( aShape );
1445     }
1446     else
1447     {
1448       TopExp_Explorer ex(aShape, TopAbs_EDGE);
1449       for(; ex.More(); ex.Next() )
1450         anEdgeMap.Add( ex.Current() );
1451     }
1452
1453     if ( anEdgeMap.Extent() == 0 )
1454       continue;
1455
1456     char newname[1024];
1457     strcpy(newname,argv[i]);
1458     char* p = newname;
1459     while (*p != '\0') p++;
1460     *p = '_';
1461     p++;
1462
1463     Standard_Integer nbEdge = 1;
1464     for(it.Initialize(anEdgeMap); it.More(); it.Next())
1465     {
1466       BRep_Tool::PolygonOnTriangulation(TopoDS::Edge(it.Key()), aPoly, aT, aLoc);
1467       if ( aT.IsNull() || aPoly.IsNull() )
1468         continue;
1469       
1470       const TColgp_Array1OfPnt&      Nodes   = aT->Nodes();
1471       const TColStd_Array1OfInteger& Indices = aPoly->Nodes();
1472       const Standard_Integer         nbnodes = Indices.Length();
1473
1474       for( Standard_Integer j = 1; j <= nbnodes; j++ )
1475       {
1476         gp_Pnt P3d = Nodes(Indices(j));
1477         if( !aLoc.IsIdentity() )
1478           P3d.Transform(aLoc.Transformation());
1479
1480         if( anEdgeMap.Extent() > 1 )
1481           Sprintf(p,"%d_%d",nbEdge,j);
1482         else
1483           Sprintf(p,"%d",j);
1484               DBRep::Set( newname, BRepBuilderAPI_MakeVertex(P3d) );
1485               di.AppendElement(newname);
1486       }
1487       nbEdge++;
1488     }
1489   }
1490   return 0;
1491 }
1492
1493 //=======================================================================
1494 void  MeshTest::Commands(Draw_Interpretor& theCommands)
1495 //=======================================================================
1496 {
1497   Draw::Commands(theCommands);
1498   BRepTest::AllCommands(theCommands);
1499   GeometryTest::AllCommands(theCommands);
1500   MeshTest::PluginCommands(theCommands);
1501   const char* g;
1502
1503   g = "Mesh Commands";
1504
1505   theCommands.Add("incmesh","incmesh shape deflection [inParallel (0/1) : 0 by default]",__FILE__, incrementalmesh, g);
1506   theCommands.Add("MemLeakTest","MemLeakTest",__FILE__, MemLeakTest, g);
1507   theCommands.Add("fastdiscret","fastdiscret shape deflection [shared [nbiter]]",__FILE__, fastdiscret, g);
1508   theCommands.Add("mesh","mesh result Shape deflection [save partage]",__FILE__, triangule, g);
1509   theCommands.Add("addshape","addshape meshname Shape [deflection]",__FILE__, addshape, g);
1510   //theCommands.Add("smooth","smooth meshname",__FILE__, smooth, g);
1511   //theCommands.Add("edges","edges mesh shape, highlight the edges",__FILE__,edges, g);
1512   theCommands.Add("vertices","vertices mesh shape, highlight the vertices",__FILE__,vertices, g);
1513   theCommands.Add("medge","medge mesh [-]index (0 to clear all)",__FILE__,medge, g);
1514   theCommands.Add("mvertex","mvertex mesh [-]index (0 to clear all)",__FILE__,mvertex, g);
1515   theCommands.Add("triangle","triangle mesh [-]index (0 to clear all)",__FILE__,triangle, g);
1516   //theCommands.Add("dumpvertex","dumpvertex mesh [index]",__FILE__,dumpvertex, g);
1517   //theCommands.Add("dumpedge","dumpedge mesh [index]",__FILE__,dumpedge, g);
1518   //theCommands.Add("dumptriangle","dumptriangle mesh [index]",__FILE__,dumptriangle, g);
1519
1520   theCommands.Add("tri2d", "tri2d facename",__FILE__, tri2d, g);
1521   theCommands.Add("trinfo","trinfo name, print triangles information on objects",__FILE__,trianglesinfo,g);
1522   theCommands.Add("veriftriangles","veriftriangles name, verif triangles",__FILE__,veriftriangles,g);
1523   theCommands.Add("wavefront","wavefront name",__FILE__, wavefront, g);
1524   theCommands.Add("onetriangulation","onetriangulation name",__FILE__, onetriangulation, g);
1525   theCommands.Add("triepoints", "triepoints shape1 [shape2 ...]",__FILE__, triedgepoints, g);
1526
1527 #if 0
1528   theCommands.Add("extrema","extrema ",__FILE__, extrema, g);
1529   theCommands.Add("vb","vb ",__FILE__, vb, g);
1530 #endif
1531 }