Allink  v0.1
ElPolyDrawCGAL.cpp
1 /***********************************************************************
2 ElPoly:This progam provide a graphical visualisation of the data
3 opend by VarData using openGL glut. The most important option are
4 the possibility of changing the backfold of the polymers with 'c',
5 see the subsequent file in the list with '>', see the bond with 'b'.
6 Copyright (C) 2008-2010 by Giovanni Marelli <sabeiro@virgilio.it>
7 
8 
9 This program is free software; you can redistribute it and/or modify
10 it under the terms of the GNU General Public License as published by
11 the Free Software Foundation; either version 2 of the License, or
12 (at your option) any later version.
13 
14 This program is distributed in the hope that it will be useful,
15 but WITHOUT ANY WARRANTY; without even the implied warranty of
16 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
17 GNU General Public License for more details.
18 
19 You should have received a copy of the GNU General Public License
20 along with this program; If not, see <http://www.gnu.org/licenses/>.
21 ***********************************************************************/
22 #include "ElPoly.h"
23 
24 #ifdef __glut_h__
25 extern Draw *Dr;
26  //#ifndef __CGAL_h__
27 #ifdef USE_CGAL
28 //#include <CGAL/basic.h>
29 //#include <CGAL/Cartesian.h>
30 #include <CGAL/Exact_predicates_inexact_constructions_kernel.h>
31 #include <CGAL/Surface_mesh_default_triangulation_3.h>
32 #include <CGAL/Complex_2_in_triangulation_3.h>
33  //#include <CGAL/Implicit_surface_3.h>
34 #include <CGAL/Delaunay_triangulation_2.h>
35 #include <CGAL/Triangulation_euclidean_traits_xy_3.h>
36 //using namespace CGAL;
37 
38 //default triangulation for Surface_mesher
39 typedef CGAL::Surface_mesh_default_triangulation_3 Tr;
40 
41 //c2t3
42 typedef CGAL::Complex_2_in_triangulation_3<Tr> C2t3;
43 
44 typedef Tr::Geom_traits GT;
45 
46 typedef GT::Sphere_3 Sphere_3;
47 //typedef GT::Point_3 Point_3;
48 typedef GT::FT FT;
49 
50 
51 typedef CGAL::Exact_predicates_inexact_constructions_kernel K;
52 typedef CGAL::Triangulation_euclidean_traits_xy_3<K> Gt;
53 typedef CGAL::Delaunay_triangulation_2<Gt> Delaunay;
54 // typedef Triangulation_default_data_structure_2<Gt,Vb,Fb > Tds;
55 // typedef Triangulation_2<Gt,Tds> Triangulation;
56 template < class Gt >
57 class My_face_base : public CGAL::Triangulation_face_base_2<Gt>
58 {
59 public:
60  CGAL::Color color;
61  My_face_base() :
62  CGAL::Triangulation_face_base_2<Gt>() {}
63  My_face_base(void* v0, void* v1, void* v2) :
64  CGAL::Triangulation_face_base_2<Gt>(v0,v1,v2) {}
65  My_face_base(void* v0, void* v1, void* v2, void* n0, void* n1, void* n2) :
66  CGAL::Triangulation_face_base_2<Gt>(v0,v1,v2,n0,n1,n2) {}
67 };
68 typedef struct {double r;double g;double b;double a;} COLOR;
69 
71  glLineWidth(.5);
72  Delaunay dtUp;
73  Delaunay dtDown;
74  //vector <COLOR> HueUp;
75  //vector <COLOR> HueDown;
76  COLOR ColChain;
77  ColChain.a = 1.;
78  double AreaMean = Gen->Edge[CLat1]*Gen->Edge[CLat2]*.5/(double)Gen->NChain;
79  Delaunay::Finite_vertices_iterator vUp = dtUp.finite_vertices_begin();
80  for(int c=0;c<Gen->NChain;c++){
81  int Chc = Ch[c].Type;
82  if(Ch[c].Pos[CNorm] > .7*Gen->Edge[CNorm]) continue;
83  if(Ch[c].Pos[CNorm] < .3*Gen->Edge[CNorm]) continue;
84  if(!CHAIN_IF_TYPE(Chc,NChType)) continue;
85  ColChain.r = 0.;
86  if(CHAIN_IF_TYPE(Chc,CHAIN_UP))
87  ColChain.r = 1.;
88  ColChain.b = 0.;
89  if(CHAIN_IF_TYPE(Chc,CHAIN_FLABBY))
90  ColChain.b = 1.;
91  ColChain.g=.7;
92  if(CHAIN_IF_TYPE(Chc,CHAIN_TILTED))
93  ColChain.g = .5;
94  if(CHAIN_IF_TYPE(Ch[c].Type,CHAIN_UP) ){
95  K::Point_3 ChPos(Ch[c].Pos[CLat1],Ch[c].Pos[CLat2],Ch[c].Pos[CNorm]);
96  dtUp.insert(ChPos);
97  //HueUp.push_back(ColChain);
98  }
99  else{
100  K::Point_3 ChPos(Ch[c].Pos[CLat1],Ch[c].Pos[CLat2],Ch[c].Pos[CNorm]);
101  dtDown.insert(ChPos);
102  //HueDown.push_back(ColChain);
103  }
104  }
105  Delaunay::Face_iterator fcTrUp = dtUp.finite_faces_begin();
106  Delaunay::Face_iterator fcTrDown = dtDown.finite_faces_begin();
107  glDeleteLists(Dr->Particles,1);
108  Dr->Particles = glGenLists(1);
109  glNewList(Dr->Particles,GL_COMPILE);
110  for(int p=0;fcTrUp != dtUp.faces_end(); ++fcTrUp,p++){
111  Delaunay::Vertex_handle vf1 = fcTrUp->vertex(0),
112  vf2 = fcTrUp->vertex(1),vf3 = fcTrUp->vertex(2);
113  K::Point_3 pf1 = vf1->point();
114  K::Point_3 pf2 = vf2->point();
115  K::Point_3 pf3 = vf3->point();
116  Vettore v1(pf1.x(),pf1.y(),pf1.z());
117  Vettore v2(pf2.x(),pf2.y(),pf2.z());
118  Vettore v3(pf3.x(),pf3.y(),pf3.z());
119  Vettore vN(3);
120  v1.Mult(InvScaleUn);
121  v2.Mult(InvScaleUn);
122  v3.Mult(InvScaleUn);
123  vN = (v1-v2) ^ (v3-v2);
124  //if(vN.Norm() > 2.*AreaMean) continue;
125  glColor4f(0.,.7,0.,1.);
126  //glColor4f(HueUp[p].r,HueUp[p].g,HueUp[p].b,HueUp[p].a);
127  DrTria(&v1,&v2,&v3,&vN);
128  glColor4f(0.,.0,0.,1.);
129  DrTriaContour(&v1,&v2,&v3);
130  }
131  for(int p=0;fcTrDown != dtDown.faces_end(); ++fcTrDown,p++){
132  Delaunay::Vertex_handle vf1 = fcTrDown->vertex(0),
133  vf2 = fcTrDown->vertex(1),vf3 = fcTrDown->vertex(2);
134  K::Point_3 pf1 = vf1->point();
135  K::Point_3 pf2 = vf2->point();
136  K::Point_3 pf3 = vf3->point();
137  Vettore v1(pf1.x(),pf1.y(),pf1.z());
138  Vettore v2(pf2.x(),pf2.y(),pf2.z());
139  Vettore v3(pf3.x(),pf3.y(),pf3.z());
140  Vettore vN(3);
141  v1.Mult(InvScaleUn);
142  v2.Mult(InvScaleUn);
143  v3.Mult(InvScaleUn);
144  vN = (v1-v2) ^ (v3-v2);
145  //if(vN.Norm() > 2.*AreaMean) continue;
146  glColor4f(0.,.7,0.,1.);
147  // glColor4f(HueDown[p].r,HueDown[p].g,HueDown[p].b,HueDown[p].a);
148  DrTria(&v1,&v2,&v3,&vN);
149  glColor4f(0.,.0,0.,1.);
150  DrTriaContour(&v1,&v2,&v3);
151  }
152  DrNano();
153  glEndList();
154 }
155 
156 //#include <CGAL/Exact_predicates_inexact_constructions_kernel.h>
157 //#include <CGAL/make_skin_surface_mesh_3.h>
158 #include <CGAL/Skin_surface_3.h>
159 #include <list>
160 #include <CGAL/mesh_skin_surface_3.h>
161 #include <CGAL/subdivide_skin_surface_mesh_3.h>
162 //#include <CGAL/Skin_surface_3.h>
163 #include <CGAL/Polyhedron_3.h>
164 #include <fstream>
165 #include "skin_surface_writer.h"
166 #include <CGAL/IO/Polyhedron_iostream.h>
167 #include <CGAL/subdivide_skin_surface_mesh_3.h>
168 #include <CGAL/Skin_surface_refinement_policy_3.h>
169 //#include <CGAL/mesh_skin_surface_3.h>
170 //#include <CGAL/subdivide_skin_surface_mesh_3.h>
171 //typedef CGAL::Exact_predicates_inexact_constructions_kernel K;
172 typedef CGAL::Skin_surface_traits_3<K> Traits;
173 typedef CGAL::Skin_surface_3<Traits> Skin_surface_3;
174 typedef Skin_surface_3::FT FT;
175 typedef Skin_surface_3::Weighted_point Weighted_point;
176 typedef Weighted_point::Point Bare_point;
177 typedef CGAL::Polyhedron_3<K,
178  CGAL::Skin_surface_polyhedral_items_3<Skin_surface_3> > Polyhedron;
179 typedef Polyhedron::Traits::Vector_3 Vector;
180 //CGAL::Skin_surface_refinement_policy_3<SkinSurface, Polyhedron> policy(skin);
181 #include <CGAL/Skin_surface_refinement_policy_3.h>
182 // typedef K::Point_3 Bare_point;
183 // typedef CGAL::Weighted_point<Bare_point,K::RT> Weighted_point;
184 // typedef CGAL::Polyhedron_3<K> Polyhedron;
185 
186 
188  list <Weighted_point> WeiPoint;
189  double shrinkfactor = 1.0;
190 
191  for(int c=0;c<4;c++){//Gen->NChain;c++){
192  Bare_point ChPoint(Ch[c].Pos[0]*InvScaleUn,Ch[c].Pos[1]*InvScaleUn,Ch[c].Pos[2]*InvScaleUn);
193  WeiPoint.push_front(Weighted_point(ChPoint, 0.5));
194  }
195  Polyhedron Polyhe;
196  //CGAL::make_skin_surface_mesh_3(Polyhe, WeiPoint.begin(), WeiPoint.end(), shrinkfactor);
197  Skin_surface_3 skin_surface(WeiPoint.begin(),WeiPoint.end(),shrinkfactor);
198  CGAL::mesh_skin_surface_3(skin_surface,Polyhe);
199  CGAL::subdivide_skin_surface_mesh_3(skin_surface,Polyhe);
200  std::ofstream out("mesh.off");
201  out << Polyhe;
202  // Polyhedron::Facet_iterator fcUp = Polyhe.facet_begin();
203  // for(;fcUp != Polyhe.faces_end(); ++fcUp){
204  // Halfedge_around_facet_circulator heUp = fcUp.halfedge();
205  glDeleteLists(Dr->Particles,1);
206  Dr->Particles = glGenLists(1);
207  glNewList(Dr->Particles,GL_COMPILE);
208  glColor4f(0.0,0.0,0.0,1.0);
209  double Pos[3];
210  // for (Polyhedron::Vertex_iterator vit = Polyhe.vertices_begin();vit != Polyhe.vertices_end(); vit ++) {
211  //Vector n = policy.normal(vit);
212  //n = n/sqrt(n*n);
213  // out << vit->point() << " " << n << std::endl;
214  // Polyhedron::Halfedge_iterator heUp = Polyhe.halfedges_begin();
215  // for(;heUp != Polyhe.halfedges_end(); ++heUp){
216  // //Polyhedron::Halfedge_handle Half = *heUp;
217  // Polyhedron::Vertex_handle veUp = heUp->vertex();
218 //K::Point_3 pf1 = vit->point();
219 
220 
221  // for(Polyhedron::Facet_iterator fi = Polyhe.facets_begin();fi != Polyhe.facets_end(); ++fi) {
222  // Skin_surface_3::HFC hc = fi->facet_begin();
223  // Polyhedron::HFC hc_end = hc;
224  // glPushMatrix();//Particle
225  // glBegin(GL_LINES);
226  // do {
227  // Polyhedron::Vertex_handle vh = (*hc).vertex();
228  // K::Point_3 pf1 = vh->point();
229  // Pos[0] = pf1.x();Pos[1] = pf1.y();Pos[2] = pf1.z();
230  // glVertex3d(Pos[0],Pos[1],Pos[2]);
231  // } while (++hc != hc_end);
232  // glEnd();
233  // glPopMatrix();//Particle
234 
235  // }
236  glEndList();
237 
238  // Tr tr; // 3D-Delaunay triangulation
239  // C2t3 c2t3 (tr); // 2D-complex in 3D-Delaunay triangulation
240 
241  // Surface_3 surface(sphere_function,Sphere_3(CGAL::ORIGIN, 2.));
242  // Surface_mesh_default_criteria_3<Tr> criteria(30., 0.1,0.1);
243  // // meshing surface
244  // make_surface_mesh(c2t3, surface, criteria, CGAL::Non_manifold_tag());
245 
246  // std::cout << "Final number of points: " << tr.number_of_vertices() << "\n";
247 
248  // DT dt;
249  // for(int c=0;c<Gen->NChain;c++){
250  // if(CHAIN_IF_TYPE(Ch[c].Type,CHAIN_UP) )continue;
251  // Point_3<K> ChPos(Ch[c].Pos[CLat1],Ch[c].Pos[CLat2],Ch[c].Pos[CNorm]);
252  // dt.insert(ChPos);
253  // }
254  // Face_iterator fcTr = dt.finite_faces_begin();
255  // glDeleteLists(Dr->Particles,1);
256  // Dr->Particles = glGenLists(1);
257  // glNewList(Dr->Particles,GL_COMPILE);
258  // for(;fcTr != dt.faces_end(); ++fcTr){
259  // Vertex_handle vf1 = fcTr->vertex(0),
260  // vf2 = fcTr->vertex(1),vf3 = fcTr->vertex(2);
261  // Point pf1 = vf1->point();
262  // Point pf2 = vf2->point();
263  // Point pf3 = vf3->point();
264  // Vettore v1(pf1.x() - pf2.x(),pf1.y() - pf2.y(),pf1.z()-pf2.z());
265  // Vettore v2(pf3.x() - pf2.x(),pf3.y() - pf2.y(),pf1.z()-pf2.z());
266  // Vettore vN(3);
267  // vN = v1 ^ v2;
268  // DrTira(v1,v2,v3,vN);
269 
270  // }
271  // glEndList();
272 }
273 #else
274 void ElPoly::DrTriangulate(){
275  printf("DefSurface without CGAL not implemented\n");
276 }
277 #endif // __CGAL_h__
278 #endif //__glut_h__
Draw provides the basic configuration of the openGL libraries used in every derived program...
Definition: Draw.h:15
Geometrical operations on vectors.
Definition: MatematicaVect.h:9
void DefineSkin(int NSample)
Find the covering surface for given points.
void DrTriangulate()
Triangulate a surface.
GLuint Particles
Refers to the list of the total position of the particles which will be generated in another program...
Definition: Draw.h:200