Allink  v0.1
ElPolyDrawTriangCGAL.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  //#include <stdafx.h>
24 
25 #ifdef __glut_h__
26 extern Draw *Dr;
27  //#ifndef __CGAL_h__
28 #ifdef USE_CGAL
29 //#include <CGAL/basic.h>
30 //#include <CGAL/Cartesian.h>
31 #include <CGAL/Exact_predicates_inexact_constructions_kernel.h>
32 #include <CGAL/Surface_mesh_default_triangulation_3.h>
33 #include <CGAL/Complex_2_in_triangulation_3.h>
34  //#include <CGAL/Implicit_surface_3.h>
35 #include <CGAL/Delaunay_triangulation_2.h>
36 #include <CGAL/Triangulation_euclidean_traits_xy_3.h>
37 //using namespace CGAL;
38 
39 //default triangulation for Surface_mesher
40 typedef CGAL::Surface_mesh_default_triangulation_3 Tr;
41 
42 //c2t3
43 typedef CGAL::Complex_2_in_triangulation_3<Tr> C2t3;
44 
45 typedef Tr::Geom_traits GT;
46 
47 typedef GT::Sphere_3 Sphere_3;
48 //typedef GT::Point_3 Point_3;
49 typedef GT::FT FT;
50 
51 
52 typedef CGAL::Exact_predicates_inexact_constructions_kernel K;
53 typedef CGAL::Triangulation_euclidean_traits_xy_3<K> Gt;
54 typedef CGAL::Delaunay_triangulation_2<Gt> Delaunay;
55 // typedef Triangulation_default_data_structure_2<Gt,Vb,Fb > Tds;
56 // typedef Triangulation_2<Gt,Tds> Triangulation;
57 template < class Gt >
58 class My_face_base : public CGAL::Triangulation_face_base_2<Gt>
59 {
60 public:
61  CGAL::Color color;
62  My_face_base() :
63  CGAL::Triangulation_face_base_2<Gt>() {}
64  My_face_base(void* v0, void* v1, void* v2) :
65  CGAL::Triangulation_face_base_2<Gt>(v0,v1,v2) {}
66  My_face_base(void* v0, void* v1, void* v2, void* n0, void* n1, void* n2) :
67  CGAL::Triangulation_face_base_2<Gt>(v0,v1,v2,n0,n1,n2) {}
68 };
69 typedef struct {double r;double g;double b;double a;} COLOR;
70 
72  if(Dr->IfMaterial){
73  glEnable(GL_LIGHTING);
74  glEnable( GL_LIGHT0 );
75  }
76  else
77  glDisable(GL_LIGHTING);
78  Delaunay dtUp;
79  Delaunay dtDown;
80  //vector <COLOR> HueUp;
81  //vector <COLOR> HueDown;
82  COLOR ColChain;
83  GLfloat Color[4];
84  ColChain.a = 1.;
85  double AreaMean = pEdge(CLat1)*pEdge(CLat2)*.5/(double)pNChain();
86  BfDefChain();
87  double Pos[3];
88  Delaunay::Finite_vertices_iterator vUp = dtUp.finite_vertices_begin();
89  for(int c=0;c<pNChain();c++){
90  int Chc = Ch[c].Type;
91  //if(Ch[c].Pos[CNorm] > .7*pEdge(CNorm)) continue;
92  //if(Ch[c].Pos[CNorm] < .3*pEdge(CNorm)) continue;
93  if(!CHAIN_IF_TYPE(Chc,NChType)) continue;
94  ColChain.r = 0.;
95  if(CHAIN_IF_TYPE(Chc,CHAIN_UP))
96  ColChain.r = 1.;
97  ColChain.b = 0.;
98  if(CHAIN_IF_TYPE(Chc,CHAIN_FLABBY))
99  ColChain.b = 1.;
100  ColChain.g=.7;
101  if(CHAIN_IF_TYPE(Chc,CHAIN_TILTED))
102  ColChain.g = .5;
103  for(int d=0;d<3;d++)
104  Pos[d] = Ch[c].Pos[d] - floor(Ch[c].Pos[d]*pInvEdge(d))*pEdge(d);
105  if(CHAIN_IF_TYPE(Ch[c].Type,CHAIN_UP) ){
106  K::Point_3 ChPos(Pos[0],Pos[1],Pos[2]);
107  dtUp.insert(ChPos);
108  //HueUp.push_back(ColChain);
109  }
110  else{
111  K::Point_3 ChPos(Pos[0],Pos[1],Pos[2]);
112  dtDown.insert(ChPos);
113  //HueDown.push_back(ColChain);
114  }
115  }
116  Delaunay::Face_iterator fcTrUp = dtUp.finite_faces_begin();
117  Delaunay::Face_iterator fcTrDown = dtDown.finite_faces_begin();
118  glDeleteLists(Dr->Particles,1);
119  Dr->Particles = glGenLists(1);
120  glNewList(Dr->Particles,GL_COMPILE);
121  for(int p=0;fcTrUp != dtUp.faces_end(); ++fcTrUp,p++){
122  Delaunay::Vertex_handle vf1 = fcTrUp->vertex(0),
123  vf2 = fcTrUp->vertex(1),vf3 = fcTrUp->vertex(2);
124  K::Point_3 pf1 = vf1->point();
125  K::Point_3 pf2 = vf2->point();
126  K::Point_3 pf3 = vf3->point();
127  Vettore v1(pf1.x(),pf1.y(),pf1.z());
128  Vettore v2(pf2.x(),pf2.y(),pf2.z());
129  Vettore v3(pf3.x(),pf3.y(),pf3.z());
130  Vettore vN(3);
131  v1.Mult(InvScaleUn);
132  v2.Mult(InvScaleUn);
133  v3.Mult(InvScaleUn);
134  vN = (v1-v2) ^ (v3-v2);
135  glEnable(GL_LIGHTING);
136  double Depth = pf1.z()*pInvEdge(CNorm)*Saturation+ExtParam;
137  Dr->DepthMap(Depth,Color);
138  glMaterialfv(GL_FRONT_AND_BACK,GL_AMBIENT,Color);
139  glMaterialfv(GL_FRONT_AND_BACK,GL_DIFFUSE,Color);
140  glColor4fv(Color);
141  //if(vN.Norm() > 2.*AreaMean) continue;
142  DrTria(&v1,&v2,&v3,&vN);
143  glDisable(GL_LIGHTING);
144  glColor4f(0.,0.,1.,1.);
145  DrTriaContour(&v1,&v2,&v3);
146  }
147  for(int p=0;fcTrDown != dtDown.faces_end(); ++fcTrDown,p++){
148  Delaunay::Vertex_handle vf1 = fcTrDown->vertex(0),
149  vf2 = fcTrDown->vertex(1),vf3 = fcTrDown->vertex(2);
150  K::Point_3 pf1 = vf1->point();
151  K::Point_3 pf2 = vf2->point();
152  K::Point_3 pf3 = vf3->point();
153  Vettore v1(pf1.x(),pf1.y(),pf1.z());
154  Vettore v2(pf2.x(),pf2.y(),pf2.z());
155  Vettore v3(pf3.x(),pf3.y(),pf3.z());
156  Vettore vN(3);
157  v1.Mult(InvScaleUn);
158  v2.Mult(InvScaleUn);
159  v3.Mult(InvScaleUn);
160  vN = (v1-v2) ^ (v3-v2);
161  glEnable(GL_LIGHTING);
162  double Depth = pf1.z()*pInvEdge(CNorm)*Saturation+ExtParam;
163  Dr->DepthMap(Depth,Color);
164  glMaterialfv(GL_FRONT_AND_BACK,GL_AMBIENT,Color);
165  glMaterialfv(GL_FRONT_AND_BACK,GL_DIFFUSE,Color);
166  glColor4fv(Color);
167  //if(vN.Norm() > 2.*AreaMean) continue;
168  // glColor4f(HueDown[p].r,HueDown[p].g,HueDown[p].b,HueDown[p].a);
169  DrTria(&v1,&v2,&v3,&vN);
170  glDisable(GL_LIGHTING);
171  glColor4f(1.,0.,0.,1.);
172  DrTriaContour(&v1,&v2,&v3);
173  }
174  glDisable(GL_LIGHTING);
175  for(int n=0;n<pNNano();n++) DrawNano(n);
176  glEndList();
177 }
178 #include <CGAL/Triangulation_3.h>
179 typedef CGAL::Triangulation_3<K> Triang;
181  Delaunay dtUp;
182  COLOR ColChain;
183  ColChain.a = 1.;
184  double AreaMean = pEdge(CLat1)*pEdge(CLat2)*.5/(double)pNChain();
185  std::list<Triang::Point> LipPos;
186  Delaunay::Finite_vertices_iterator vUp = dtUp.finite_vertices_begin();
187  for(int b=0,cOff=0;b<pNBlock();b++,cOff+=Block[b].NChain){
188  for(int c=cOff;c<pNChain(b)+cOff;c+=20){
189  int Chc = Ch[c].Type;
190  if(!CHAIN_IF_TYPE(Chc,NChType)) continue;
191  ColChain.r = 0.;
192  if(CHAIN_IF_TYPE(Chc,CHAIN_UP))
193  ColChain.r = 1.;
194  ColChain.b = 0.;
195  if(CHAIN_IF_TYPE(Chc,CHAIN_FLABBY))
196  ColChain.b = 1.;
197  ColChain.g=.7;
198  if(CHAIN_IF_TYPE(Chc,CHAIN_TILTED))
199  ColChain.g = .5;
200  Triang::Point LipPoint(Ch[c].Pos[CLat1],Ch[c].Pos[CLat2],Ch[c].Pos[CNorm]);
201  LipPos.push_front(LipPoint);
202  //HueUp.push_back(ColChain);
203  }
204  }
205  Triang Tr(LipPos.begin(),LipPos.end());
206  Triang::size_type NVert = Tr.number_of_vertices();
207  std::vector<Triang::Point> V(3);
208  assert(Tr.is_valid());
209  // Triang::Locate_type lt;
210  // int li,lj;
211  // Triang::Point p(0,0,0);
212  // Triang::Cell_handle c = Tr.locate(p,lt,li,lj);
213  // assert( lt == Triang::VERTEX);
214  // assert(c->vertex(li)->point() == p);
215  // Triang::Vertex_handle v = c->vertex( (li+1)&3);
216  // Triang::Cell_handle nc = c->neighbor(li);
217  // int nli;
218  // assert(nc->has_vertex(v,nli));
219  glDeleteLists(Dr->Particles,1);
220  Dr->Particles = glGenLists(1);
221  glNewList(Dr->Particles,GL_COMPILE);
222  //Triang::Cell_iterator CIt = Tr.finite_cells_begin();
223  Triang::Finite_cells_iterator CIt = Tr.finite_cells_begin();
224  for(int p=0;CIt != Tr.finite_cells_end();++CIt,p++){
225  Triang::Vertex_handle vf1 = CIt->vertex(0),vf2 = CIt->vertex(1),vf3 = CIt->vertex(2);
226  K::Point_3 pf1 = vf1->point();
227  K::Point_3 pf2 = vf2->point();
228  K::Point_3 pf3 = vf3->point();
229  Vettore v1(pf1.x(),pf1.y(),pf1.z());
230  Vettore v2(pf2.x(),pf2.y(),pf2.z());
231  Vettore v3(pf3.x(),pf3.y(),pf3.z());
232  Vettore vN(3);
233  v1.Mult(InvScaleUn);
234  v2.Mult(InvScaleUn);
235  v3.Mult(InvScaleUn);
236  vN = (v1-v2) ^ (v3-v2);
237  //if(vN.Norm() > 2.*AreaMean) continue;
238  double Sfumatura = .3*Mat->Casuale();
239  glColor4f(0.1,.4+Sfumatura,0.2,1.);
240  //glColor4f(HueUp[p].r,HueUp[p].g,HueUp[p].b,HueUp[p].a);
241  DrTria(&v1,&v2,&v3,&vN);
242  glColor4f(0.,.0,0.,1.);
243  DrTriaContour(&v1,&v2,&v3);
244  }
245  for(int n=0;n<pNNano();n++) DrawNano(n);
246  glEndList();
247 }
248 #include <CGAL/Surface_mesh_default_triangulation_3.h>
249 #include <CGAL/Complex_2_in_triangulation_3.h>
250 #include <CGAL/make_surface_mesh.h>
251 #include <CGAL/Implicit_surface_3.h>
252 
253 // default triangulation for Surface_mesher
254 typedef CGAL::Surface_mesh_default_triangulation_3 Tr;
255 
256 // c2t3
257 typedef CGAL::Complex_2_in_triangulation_3<Tr> C2t3;
258 
259 typedef Tr::Geom_traits GT;
260 typedef GT::Sphere_3 Sphere_3;
261 typedef GT::Point_3 Point_3;
262 typedef GT::FT FT;
263 
264 typedef FT (*Function)(Point_3);
265 
266 typedef CGAL::Implicit_surface_3<GT, Function> Surface_3;
267 
268 FT sphere_function (Point_3 p) {
269  const FT x2=p.x()*p.x(), y2=p.y()*p.y(), z2=p.z()*p.z();
270  return x2+y2+z2-1;
271 }
273  Tr tr; // 3D-Delaunay triangulation
274  C2t3 c2t3 (tr); // 2D-complex in 3D-Delaunay triangulation
275  Surface_3 surface(sphere_function,Sphere_3(CGAL::ORIGIN, 2.));
276  CGAL::Surface_mesh_default_criteria_3<Tr> criteria(30., // angular bound
277  0.1, // radius bound
278  0.1); // distance bound
279  CGAL::make_surface_mesh(c2t3, surface, criteria, CGAL::Non_manifold_tag());
280  glDeleteLists(Dr->Particles,1);
281  Dr->Particles = glGenLists(1);
282  glNewList(Dr->Particles,GL_COMPILE);
283  Tr::Finite_cells_iterator CIt = tr.finite_cells_begin();
284  for(int p=0;CIt != tr.finite_cells_end();++CIt,p++){
285  Tr::Vertex_handle vf1 = CIt->vertex(0),vf2 = CIt->vertex(1),vf3 = CIt->vertex(2);
286  K::Point_3 pf1 = vf1->point();
287  K::Point_3 pf2 = vf2->point();
288  K::Point_3 pf3 = vf3->point();
289  Vettore v1(pf1.x(),pf1.y(),pf1.z());
290  Vettore v2(pf2.x(),pf2.y(),pf2.z());
291  Vettore v3(pf3.x(),pf3.y(),pf3.z());
292  Vettore vN(3);
293  v1.Mult(InvScaleUn);
294  v2.Mult(InvScaleUn);
295  v3.Mult(InvScaleUn);
296  vN = (v1-v2) ^ (v3-v2);
297  //if(vN.Norm() > 2.*AreaMean) continue;
298  double Sfumatura = .3*Mat->Casuale();
299  glColor4f(0.1,.4+Sfumatura,0.2,1.);
300  //glColor4f(HueUp[p].r,HueUp[p].g,HueUp[p].b,HueUp[p].a);
301  DrTria(&v1,&v2,&v3,&vN);
302  glColor4f(1.,.0,0.,1.);
303  DrTriaContour(&v1,&v2,&v3);
304  }
305  for(int n=0;n<pNNano();n++) DrawNano(n);
306  glEndList();
307 }
308 #else
309 void ElPoly::DrTriangulate(){
310  printf("DefSurface without CGAL not implemented\n");
311 }
312 void ElPoly::DrMesh(){
313  printf("DefSurface without CGAL not implemented\n");
314 }
315 void ElPoly::DrCells(){
316  printf("DefSurface without CGAL not implemented\n");
317 }
318 #endif // __CGAL_h__
319 #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 DrGenMesh()
Generate a mesh from a function.
void DepthMap(double Val, GLfloat *Color)
Pointer to a generic function.
Definition: Draw.h:113
void DrTriangulate()
Triangulate a surface.
void DrCells()
Construct cells from the lipid positions.
int IfMaterial
Activate the illumination for a specific material.
Definition: Draw.h:240
void DrMesh()
Build a mesh from the lipid positions.
GLuint Particles
Refers to the list of the total position of the particles which will be generated in another program...
Definition: Draw.h:200