Allink  v0.1
ElPolyDrawSurf.cpp
1 /*
2 * ElPolyDrawSurf - Drawing of the surfaces defined by the chain position
3  opyright (C) 2008-2010 by Giovanni Marelli <sabeiro@virgilio.it>
4 
5 
6  This program is free software; you can redistribute it and/or modify
7  it under the terms of the GNU General Public License as published by
8  the Free Software Foundation; either version 2 of the License, or
9  (at your option) any later version.
10 
11  This program is distributed in the hope that it will be useful,
12  but WITHOUT ANY WARRANTY; without even the implied warranty of
13  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14  GNU General Public License for more details.
15 
16  You should have received a copy of the GNU General Public License
17  along with this program; If not, see <http://www.gnu.org/licenses/>.
18 ***********************************************************************/
19 #include "ElPoly.h"
20 #include <vector>
21 
22 #ifdef __glut_h__
23 extern Draw *Dr;
25  glEnable(GL_LIGHTING);
26  glEnable( GL_LIGHT0 );
27  glDeleteLists(Dr->Particles,1);
28  Dr->Particles = glGenLists(1);
29  glNewList(Dr->Particles,GL_COMPILE);
30  for(int n=0;n<pNNano();n++) DrawNano(n);
31  double Cm[3];
32  for(int p=0;p<pNPart();p+=Ln[p].NLink+1){
33  //if(Pm[p].Pos[2] > .6*pEdge(2) || Pm[p].Pos[2] < .3*pEdge(2)) continue;
34  double Red = 0.0;
35  double Green = .4 + .3*Mat->Casuale();
36  double Blue = 0.0;
37  GLfloat DrMatColor [] = {Red,Green,Blue,1.};
38  glMaterialfv(GL_FRONT_AND_BACK,GL_AMBIENT,DrMatColor);
39  glMaterialfv(GL_FRONT_AND_BACK,GL_DIFFUSE,DrMatColor);
40  glColor4f(Red,Green,Blue,1.);
41  glPushMatrix();
42  glBegin(GL_POLYGON);
43  //glBegin(GL_TRIANGLE);
44  //glBegin(GL_TRIANGLE_STRIP);
45  glNormal3d(pVel(p,0)*InvScaleUn,pVel(p,1)*InvScaleUn,pVel(p,2)*InvScaleUn);
46  glVertex3d(pPos(p,0)*InvScaleUn,pPos(p,1)*InvScaleUn,pPos(p,2)*InvScaleUn*ScaleFact);
47  for(int d=0;d<3;d++) Cm[d] = pPos(p,d);
48  for(int l=0;l<Ln[p].NLink;l++){
49  int l1 = Ln[p].Link[l];
50  glNormal3d(pVel(l1,0),pVel(l1,1),pVel(l1,2));
51  glVertex3d(pPos(l1,0)*InvScaleUn,pPos(l1,1)*InvScaleUn,pPos(l1,2)*InvScaleUn*ScaleFact);
52  for(int d=0;d<3;d++){
53  Cm[d] += pPos(p,d);
54  }
55  }
56  glEnd();
57  glPopMatrix();
58  if(IfLine){
59  for(int d=0;d<3;d++){
60  Cm[d] /= (double)(Ln[p].NLink+1);
61  }
62  glColor4f(.1,.3,1.,1.);
63  glBegin(GL_LINES);
64  glNormal3d(0.,0.,1.);
65  glVertex3d(Cm[0]*InvScaleUn,Cm[1]*InvScaleUn,Cm[2]*InvScaleUn*ScaleFact);
66  glVertex3d((Cm[0]-pVel(p,0))*InvScaleUn,(Cm[1]-pVel(p,1))*InvScaleUn,(Cm[2]-pVel(p,2))*InvScaleUn*ScaleFact);
67  glEnd();
68  }
69  }
70  glEndList();
71 }
72 void ElPoly::DrNormalPoint(int p,int NEdge){
73  int link[4];
74  link[0] = p + NEdge;
75  if(link[0] >= pNPart()) link[0] -= pNPart();
76  link[1] = p + 1;
77  if(link[1] >= pNPart()) link[1] -= pNPart();
78  link[2] = p - NEdge;
79  if(link[2] < 0 ) link[2] += pNPart();
80  link[3] = p - 1;
81  if(link[3] < 0 ) link[3] += pNPart();
82  double Normals[5][3];
83  for(int l=0;l<4;l++){
84  for(int d=0;d<3;d++){
85  int d1 = (d+1)%3;
86  int d2 = (d+2)%3;
87  Normals[l][d] = (pPos(p,d1)-pPos(link[l],d2))*(pPos(p,d2)-pPos(link[l],d));
88  }
89  }
90  double Norm;
91  for(int d=0;d<3;d++){
92  Normals[4][d] = Normals[0][d] + Normals[1][d] + Normals[2][d] + Normals[3][d];
93  Norm += SQR(Normals[4][d]);
94  }
95  Norm = 1./sqrt(Norm);
96  for(int d=0;d<3;d++){
97  Normals[4][d] /= Norm;
98  }
99  glNormal3d(Normals[4][0],Normals[4][1],Normals[4][2]);
100 }
102  if(Dr->IfMaterial){
103  glEnable(GL_LIGHTING);
104  glEnable( GL_LIGHT0 );
105  }
106  else
107  glDisable(GL_LIGHTING);
108  glDeleteLists(Dr->Particles,1);
109  Dr->Particles = glGenLists(1);
110  glNewList(Dr->Particles,GL_COMPILE);
111  int NEdge = (int)sqrt(pNPart());
112  if(NEdge*NEdge != pNPart()) NEdge += 1;
113  GLfloat Color[4];
114  for(int p=0;p<pNPart();p++){
115  double Depth = pPos(p,2)*pInvEdge(CNorm)*Saturation+ExtParam;
116  Dr->DepthMap(Depth,Color);
117  glMaterialfv(GL_FRONT_AND_BACK,GL_AMBIENT,Color);
118  glMaterialfv(GL_FRONT_AND_BACK,GL_DIFFUSE,Color);
119  glColor4fv(Color);
120  int l1 = p+1;
121  if((l1%NEdge) == 0)continue;
122  int l2 = p + NEdge + 1;
123  if(l2 >= pNPart()-1)continue;
124  int l3 = p + NEdge;
125  if(l3 >= pNPart()-1)continue;
126  //glPushMatrix();
127  glBegin(GL_POLYGON);
128  DrNormalPoint(p,NEdge);
129  glVertex3d(Pm[p].Pos[0]*InvScaleUn,Pm[p].Pos[1]*InvScaleUn,Pm[p].Pos[2]*InvScaleUn*ScaleFact);
130  DrNormalPoint(l1,NEdge);
131  glVertex3d(Pm[l1].Pos[0]*InvScaleUn,Pm[l1].Pos[1]*InvScaleUn,Pm[l1].Pos[2]*InvScaleUn*ScaleFact);
132  DrNormalPoint(l2,NEdge);
133  glVertex3d(Pm[l2].Pos[0]*InvScaleUn,Pm[l2].Pos[1]*InvScaleUn,Pm[l2].Pos[2]*InvScaleUn*ScaleFact);
134  DrNormalPoint(l3,NEdge);
135  glVertex3d(Pm[l3].Pos[0]*InvScaleUn,Pm[l3].Pos[1]*InvScaleUn,Pm[l3].Pos[2]*InvScaleUn*ScaleFact);
136  glEnd();
137  //glPopMatrix();
138  }
139  glEndList();
140  glDisable(GL_LIGHTING);
141 }
142 #include "ElPolyDrawSurf.h"
143 void ElPoly::DrField(int NGrid,double IsoLevel,int nNano){
144  glEnable(GL_LIGHTING);
145  glEnable( GL_LIGHT0 );
146  glMaterialfv(GL_BACK, GL_AMBIENT, DrAmbientRed);
147  glMaterialfv(GL_BACK, GL_DIFFUSE, DrDiffuseRed);
148  glMaterialfv(GL_FRONT, GL_AMBIENT, DrAmbientBlue);
149  glMaterialfv(GL_FRONT, GL_DIFFUSE, DrDiffuseBlue);
150  glMaterialfv(GL_FRONT, GL_SPECULAR, DrSpecularWhite);
151  glMaterialf( GL_FRONT, GL_SHININESS, 25.0);
152  double CubeDist[8];
153  double EdgeVertex[12][3];
154  double EdgeNormal[12][3];
155  double InvNGrid = 1./(double)NGrid;
156  double Pos[3];
157  double Pos1[3];
158  glBegin(GL_TRIANGLES);
159  //glBegin(GL_LINES);
160  //glBegin(GL_POINTS);
161  for(int gx=0;gx<NGrid;gx++){
162  Pos[0] = gx*InvNGrid*pEdge(0);
163  for(int gy=0;gy<NGrid;gy++){
164  Pos[1] = gy*InvNGrid*pEdge(1);
165  for(int gz=0;gz<NGrid;gz++){
166  Pos[2] = gz*InvNGrid*pEdge(2);
167  for(int v=0;v<8;v++){
168  for(int d=0;d<3;d++){
169  Pos1[d] = Pos[d] + VertCube[v][d]*InvNGrid*pEdge(d);
170  }
171  CubeDist[v] = NanoDist2(Pos1,nNano);
172  }
173  int Flag = 0;
174  for(int v=0;v<8;v++){
175  if(CubeDist[v] <= IsoLevel)
176  Flag |= 1<<v;
177  }
178  int CubeShape = CubeTop[Flag];
179  if(CubeShape==0) continue;
180  for(int e=0;e<12;e++){
181  if(CubeShape & (1<<e)){
182  double Delta = CubeDist[EdgeConn[e][1]] - CubeDist[EdgeConn[e][0]];
183  double OffSet = (IsoLevel-CubeDist[EdgeConn[e][0]])/Delta;
184  if(Delta == 0.0){
185  OffSet = .5;
186  }
187  EdgeNormal[e][0] = NanoDist2(Pos[0]-0.01,Pos[1],Pos[2],nNano)
188  - NanoDist2(Pos[0]+0.01,Pos[1],Pos[2],nNano);
189  EdgeNormal[e][1] = NanoDist2(Pos[0],Pos[1]-0.01,Pos[2],nNano)
190  - NanoDist2(Pos[0],Pos[1]+0.01,Pos[2],nNano);
191  EdgeNormal[e][2] = NanoDist2(Pos[0],Pos[1],Pos[2]-0.01,nNano)
192  - NanoDist2(Pos[0],Pos[1],Pos[2]+0.01,nNano);
193  double Norm = sqrt(SQR(EdgeNormal[e][0])+SQR(EdgeNormal[e][1])+SQR(EdgeNormal[e][2]));
194  for(int d=0;d<3;d++){
195  EdgeVertex[e][d] = Pos[d] + (VertCube[EdgeConn[e][0]][d]+OffSet*EdgeDir[e][d])*InvNGrid*pEdge(d);
196  EdgeVertex[e][d] *= InvScaleUn;
197  EdgeNormal[e][d] *= 1./Norm;
198  }
199  }
200  }
201  for(int t=0;t<5;t++){
202  if(TrConnTable[Flag][3*t] < 0.) break;
203  for(int d=0;d<3;d++){
204  int v = TrConnTable[Flag][3*t+d];
205  //glColor4f(1.0,0.1,0.1,1.0);
206  glNormal3d(EdgeNormal[v][0],EdgeNormal[v][1],EdgeNormal[v][2]);
207  glVertex3d(EdgeVertex[v][0],EdgeVertex[v][1],EdgeVertex[v][2]);
208  }
209  }
210  }
211  }
212  }
213  glEnd();
214 }
216  glEnable(GL_LIGHTING);
217  glEnable( GL_LIGHT0 );
218  glDeleteLists(Dr->Particles,1);
219  Dr->Particles = glGenLists(1);
220  glNewList(Dr->Particles,GL_COMPILE);
221  //glEnable(GL_POLYGON_STIPPLE);
222  double RedMax = fabs(log(fabs(pVel(0,0))*Saturation));
223  for(int p=0;p<pNPart();p++){
224  if(RedMax < fabs(log(fabs(pVel(0,0))*Saturation)) )
225  RedMax = fabs(log(fabs(pVel(0,0))*Saturation));
226  }
227  int NDraw = 0;
228  for(int p=0;p<pNPart();p++){
229  if(pPos(p,2) < .7*pEdge(2)) continue;
230  //double Red = fabs(log(fabs(pPos(p,0))*Sat))/RedMax;
231  double Red = fabs(pVel(p,0)*Saturation/pVelMax(0));
232  double Green = pVel(p,1)/pVelMax(1)*Saturation;
233  double Blue = pVel(p,2)/pVelMax(2)*Saturation;
234  GLfloat DrMatColor [] = {Red,Green,Blue,1.};
235  glMaterialfv(GL_FRONT_AND_BACK,GL_AMBIENT,DrMatColor);
236  glMaterialfv(GL_FRONT_AND_BACK,GL_DIFFUSE,DrMatColor);
237  //if(Red<0.1) continue;
238  glColor4f(Red,Green,Blue,1.);
239  glPushMatrix();
240  glTranslated(pPos(p,0)*InvScaleUn,pPos(p,1)*InvScaleUn,pPos(p,2)*InvScaleUn*ScaleFact);
241  glCallList(Quad);
242  glPopMatrix();
243  NDraw++;
244  }
245  printf("%d\n",NDraw);
246  glEndList();
247 }
249  glEnable(GL_LIGHTING);
250  glEnable( GL_LIGHT0 );
251  //CompileList();
252  glDeleteLists(Dr->Particles,1);
253  Dr->Particles = glGenLists(1);
254  glNewList(Dr->Particles,GL_COMPILE);
255  int Incr = 1;//4;
256  double Offset = 1./(double)NEdge;
257  for(int p=0;p<pNPart();p+=Incr){
258  double Red = pVel(p,0)*Saturation;// + pPos(p,2)*Sat;
259  double Green = pVel(p,1)*Saturation;
260  double Blue = pVel(p,2)*Saturation;
261  GLfloat DrMatColor [] = {Red,Green,Blue,1.};
262  glMaterialfv(GL_FRONT_AND_BACK,GL_AMBIENT,DrMatColor);
263  glMaterialfv(GL_FRONT_AND_BACK,GL_DIFFUSE,DrMatColor);
264  //if(Red < 0.0001 && Green < 0.2 && Blue < 0.2) continue;
265  glColor4f(Red,Green,Blue,1.0);
266  glPushMatrix();
267  glTranslated(pPos(p,0)*InvScaleUn,pPos(p,1)*InvScaleUn,0.);
268  glCallList(Quad);
269  glPopMatrix();
270  // glPushMatrix();
271  // glTranslated(-pPos(p,0)*InvScaleUn-Offset,pPos(p,1)*InvScaleUn,0.);
272  // glCallList(Quad);
273  // glPopMatrix();
274  }
275  glEndList();return;
276 }
278  glEnable(GL_LIGHTING);
279  glEnable( GL_LIGHT0 );
280  glDeleteLists(Dr->Particles,1);
281  Dr->Particles = glGenLists(1);
282  glNewList(Dr->Particles,GL_COMPILE);
283  int Nx = 0;
284  int Ny = 0;
285  for(int p=0;p<pNPart();p++){
286  if(pPos(p,CLat2) < pPos(p+1,CLat2))
287  Ny++;
288  else break;
289  }
290  Ny++;
291  Nx = (int)(pNPart()/(double)Ny);
292  const int NIso = 20;
293  double Min = pPos(0,CNorm);
294  double Max = pPos(0,CNorm);
295  for(int p=0;p<pNPart();p++){
296  if(Min > pPos(p,CNorm)) Min = pPos(p,CNorm);
297  if(Max < pPos(p,CNorm)) Max = pPos(p,CNorm);
298  }
299  Min = 1./MIN(fabs(Min),1.);
300  Max = 1./MIN(1.,Max);
301  for(int p=0;p<pNPart();p++){
302  int l1 = p+1;
303  if( !(l1%Ny) ) continue;
304  int l2 = p + Ny;
305  if(l2 >= pNPart() - 1) continue;
306  int l3 = p + Ny + 1;
307  double Red = 1.;
308  double Blue = 1.;
309  double Green = 1.;
310  double Alpha = 1.;
311  if(pPos(p,CNorm) < 0.){
312  Red = 0.;
313  Green = fabs(pPos(p,CNorm)*Min)*Saturation;
314  //Green = fabs(pPos(p,2)*1.)*Sat;
315  Blue = 0.;
316  Alpha = MAX(Green-.2,0.);
317  }
318  else if(pPos(p,CNorm) > 0.){
319  Red = 0.;
320  Green = 0.;
321  Blue = (pPos(p,CNorm)*pInvEdge(CNorm))*Saturation;
322  //Blue = pPos(p,2)*1.*Sat;
323  Alpha = MAX(Blue,0.);
324  }
325  else
326  Alpha = 0.;
327  GLfloat DrMatColor [] = {Red,Green,Blue,Alpha};
328  glMaterialfv(GL_FRONT_AND_BACK,GL_AMBIENT,DrMatColor);
329  glMaterialfv(GL_FRONT_AND_BACK,GL_DIFFUSE,DrMatColor);
330  glColor4d(Red,Green,Blue,Alpha);
331  // glPushMatrix();//Particle
332  glBegin(GL_POINTS);
333  glVertex3f(pPos(p,0)*InvScaleUn,pPos(p,1)*InvScaleUn,pPos(p,2)*InvScaleUn*ScaleFact);
334  glEnd();
335  // glPopMatrix();//Particle
336  glBegin(GL_POLYGON);
337  glNormal3f(0.,0.,1.);
338  glVertex3d(pPos(p,0)*InvScaleUn,pPos(p,1)*InvScaleUn,pPos(p,2)*InvScaleUn);
339  glVertex3d(pPos(l1,0)*InvScaleUn,pPos(l1,1)*InvScaleUn,pPos(l1,2)*InvScaleUn);
340  glVertex3d(pPos(l3,0)*InvScaleUn,pPos(l3,1)*InvScaleUn,pPos(l3,2)*InvScaleUn);
341  glVertex3d(pPos(l2,0)*InvScaleUn,pPos(l2,1)*InvScaleUn,pPos(l2,2)*InvScaleUn);
342  glEnd();
343  }
344  glEndList();
345  // int Values = NEdge/4;
346  // double *Border = (double *)calloc(2*Values,sizeof(double));
347  // double *Norm = (double *)calloc(2*Values,sizeof(double));
348  // for(int p=0;p<pNPart();p++){
349  // int v = (int)(pPos(p,0)/pEdge(0)*Values);
350  // if( v >= Values) continue;
351  // if(pPos(p,1) > .5*pEdge(1)){
352  // Border[2*v] += pPos(p,1)*pPos(p,2);
353  // Norm[2*v] += pPos(p,2);
354  // }
355  // else{
356  // Border[2*v+1] += pPos(p,1)*pPos(p,2);
357  // Norm[2*v+1] += pPos(p,2];
358  // }
359  // }
360  // glColor4f(.4,0.0,.6,1.0);
361  // glBegin(GL_LINE_STRIP);
362  // glNormal3f(0.,0.,1.);
363  // for(int v=0,vv=0;v<2*Values-1;v+=2,vv++){
364  // if(!(Norm[v] > 0. ) )continue;
365  // glVertex3f(vv/(double)Values*pEdge(0)*InvScaleUn,Border[v]*InvScaleUn/Norm[v],0.);
366  // }
367  // glEnd();
368  // glBegin(GL_LINE_STRIP);
369  // for(int v=1,vv=0;v<2*Values-1;v+=2,vv++){
370  // if(!(Norm[v] > 0. ) )continue;
371  // glVertex3f(vv/(double)Values*pEdge(0)*InvScaleUn,Border[v]*InvScaleUn/Norm[v],0.);
372  // }
373  // glEnd();
374  // glBegin(GL_LINE_STRIP);
375  // for(int v=0,vv=0;v<2*Values-1;v+=2,vv++){
376  // if(!(Norm[v] > 0. ) )continue;
377  // glVertex3f(-vv/(double)Values*pEdge(0)*InvScaleUn,Border[v]*InvScaleUn/Norm[v],0.);
378  // }
379  // glEnd();
380  // glBegin(GL_LINE_STRIP);
381  // for(int v=1,vv=0;v<2*Values-1;v+=2,vv++){
382  // if(!(Norm[v] > 0. ) )continue;
383  // glVertex3f(-vv/(double)Values*pEdge(0)*InvScaleUn,Border[v]*InvScaleUn/Norm[v],0.);
384  // }
385  // glEnd();
386  // glEndList();
387 }
389  if(Dr->IfMaterial){
390  glEnable(GL_LIGHTING);
391  glEnable( GL_LIGHT0 );
392  }
393  else
394  glDisable(GL_LIGHTING);
395  glDeleteLists(Dr->Particles,1);
396  Dr->Particles = glGenLists(1);
397  glNewList(Dr->Particles,GL_COMPILE);
398  double Offset = 1./(double)NEdge;
399  for(int p=0,c=0,t=0;p<pNPart();p++){
400  if( Ln[p].NLink != 3 ) continue;
401  if(c == pChain(p)) continue;
402  if(NPType != BEAD_EVERY && pType(p) != NPType) continue;
403  c = pChain(p);
404  double Red = pVel(p,0)*Saturation/pVelMax(0);// + pPos(p,2)*Sat;
405  double Green = pVel(p,1)*Saturation/pVelMax(1);
406  double Blue = pVel(p,2)*Saturation/pVelMax(2);
407  double Alpha = 1.0;//7.+(Green+Blue);
408  if(Blue + Red + Green < .08){Blue = 1.; Red = 1.;Green = 1.;Alpha=0.;}
409  GLfloat DrMatColor [] = {Red,Green,Blue,Alpha};
410  glMaterialfv(GL_FRONT_AND_BACK,GL_AMBIENT,DrMatColor);
411  glMaterialfv(GL_FRONT_AND_BACK,GL_DIFFUSE,DrMatColor);
412  int l1 = Ln[p].Link[0];
413  //if(l1 > Gen->NPart) continue;
414  int l2 = Ln[p].Link[2];
415  //if(l2 > Gen->NPart) continue;
416  int l3 = Ln[p].Link[1];
417  //if(Green > 0.) continue;
418  // Green = 0.; Blue = 0.;
419  glColor4f(Red,Green,Blue,Alpha);
420  // if(pPos(p,2)>12.) pPos(p,2) = 12.;
421  // if(pPos(l1,2)>12.) pPos(l1,2) = 12.;
422  // if(pPos(l2,2)>12.) pPos(l2,2) = 12.;
423  // if(pPos(l3,2)>12.) pPos(l3,2) = 12.;
424  Vettore v1(pPos(p,0),pPos(p,1),pPos(p,2)*ScaleFact);
425  Vettore v2(pPos(l1,0),pPos(l1,1),pPos(l1,2)*ScaleFact);
426  Vettore v3(pPos(l2,0),pPos(l2,1),pPos(l2,2)*ScaleFact);
427  Vettore v4(pPos(l3,0),pPos(l3,1),pPos(l3,2)*ScaleFact);
428  v1.Mult(InvScaleUn);
429  v2.Mult(InvScaleUn);
430  v3.Mult(InvScaleUn);
431  v4.Mult(InvScaleUn);
432  Vettore u(3);
433  Vettore v(3);
434  Vettore vN(3);
435  vN.NormalSurf(&v1,&v2,&v3);
436  glPushMatrix();//Particle
437  glBegin(GL_QUADS);
438  glNormal3f(vN.x[0],vN.x[1],vN.x[2]);
439  u = v2-v1;v = v3-v1;vN = u ^ v;//vN.Normalize();
440  glNormal3f(vN.x[0],vN.x[1],vN.x[2]);
441  glVertex3d(v4.x[0],v4.x[1],v4.x[2]);
442  u = v1-v2;v = v3-v2;vN = u ^ v;//vN.Normalize();
443  glNormal3f(vN.x[0],vN.x[1],vN.x[2]);
444  glVertex3d(v3.x[0],v3.x[1],v3.x[2]);
445  u = v2-v3;v = v4-v3;vN = u ^ v;//vN.Normalize();
446  glNormal3f(vN.x[0],vN.x[1],vN.x[2]);
447  glVertex3d(v2.x[0],v2.x[1],v2.x[2]);
448  u = v1-v4;v = v3-v4;vN = u ^ v;//vN.Normalize();
449  glNormal3f(vN.x[0],vN.x[1],vN.x[2]);
450  glVertex3d(v1.x[0],v1.x[1],v1.x[2]);
451  glEnd();
452  glPopMatrix();//Particle
453  //continue;
454  glPushMatrix();//Particle
455  glBegin(GL_QUADS);
456  glNormal3f(vN.x[0],vN.x[1],vN.x[2]);
457  glVertex3d(-v1.x[0],v1.x[1],v1.x[2]);
458  glVertex3d(-v2.x[0],v2.x[1],v2.x[2]);
459  glVertex3d(-v3.x[0],v3.x[1],v3.x[2]);
460  glVertex3d(-v4.x[0],v4.x[1],v4.x[2]);
461  glEnd();
462  glPopMatrix();//Particle
463  }
464  glEndList();
465 }
467  glEnable(GL_LIGHTING);
468  glEnable( GL_LIGHT0 );
469  int NPoint = 32;
470  GLfloat Deltax = 1./(GLfloat) NPoint;
471  PART *PSample = (PART *)calloc(NPoint*NPoint,sizeof(PART));
472  int NSample = (int)sqrt(pNPart());
473  glColor4f(1.0,.0,1.,1.);
474  glEnableClientState(GL_VERTEX_ARRAY);
475  glEnableClientState(GL_NORMAL_ARRAY);
476  glEnable(GL_VERTEX_ARRAY);
477  glPolygonMode(GL_FRONT_AND_BACK, GL_FILL);
478  glFrontFace(GL_CCW); //Tell OGL which orientation shall be the front face
479  glDeleteLists(Dr->Particles,1);
480  Dr->Particles = glGenLists(1);
481  glNewList(Dr->Particles,GL_COMPILE);
482  //InterSpline3(Pm,PSample,NSample*NSample,NPoint*NPoint);
483  glVertexPointer(3,GL_DOUBLE,sizeof(PART),&PSample[0].Pos[0]);
484  glNormalPointer(GL_DOUBLE,sizeof(PART),&PSample[0].Vel[0]);
485  vector <GLuint> VecIndices;
486  for(int i=0;i<NPoint-1;i++){
487  for(int j=0;j<NPoint-1;j++){
488  VecIndices.push_back(i*NPoint + j);
489  VecIndices.push_back((i+1)*NPoint + j);
490  VecIndices.push_back((i+1)*NPoint + (j+1));
491 
492  VecIndices.push_back(i*NPoint + j);
493  VecIndices.push_back((i+1)*NPoint + (j+1));
494  VecIndices.push_back((i)*NPoint + (j+1));
495  }
496  }
497  int NIndex = VecIndices.size();
498  GLuint *Indices = (GLuint *)calloc(NIndex,sizeof(GLuint));
499  for(int i=0;i<NIndex;i++){
500  Indices[i] = VecIndices[i];
501  }
502  glDrawElements(GL_TRIANGLES,NIndex,GL_UNSIGNED_INT,Indices);
503  VecIndices.clear();
504  free(Indices);
505  glEndList();
506  return;
507 }
509  glPushMatrix();//Particle
510  glBegin(GL_TRIANGLES);
511  // vN->NormalSurf(v00,v01,v11);
512  // vN->Normalize();
513  // glNormal3f(vN->x[0],vN->x[1],vN->x[2]);
514  // glVertex3d(v00->x[0],v00->x[1],v00->x[2]);
515  // vN->NormalSurf(v01,v00,v11);
516  // vN->Normalize();
517  // glNormal3f(vN->x[0],vN->x[1],vN->x[2]);
518  // glVertex3d(v01->x[0],v01->x[1],v01->x[2]);
519  // vN->NormalSurf(v11,v01,v10);
520  // vN->Normalize();
521  // glNormal3f(vN->x[0],vN->x[1],vN->x[2]);
522  // glVertex3d(v11->x[0],v11->x[1],v11->x[2]);
523  // glEnd();
524  vN->NormalSurf(v00,v01,v11);
525  vN->Normalize();
526  glNormal3f(vN->x[0],vN->x[1],vN->x[2]);
527  glVertex3d(v00->x[0],v00->x[1],v00->x[2]);
528  glVertex3d(v01->x[0],v01->x[1],v01->x[2]);
529  glVertex3d(v11->x[0],v11->x[1],v11->x[2]);
530  glEnd();
531  glPopMatrix();//Particle
532  vN->NormalSurf(v00,v11,v10);
533  vN->Normalize();
534  glPushMatrix();//Particle
535  glBegin(GL_TRIANGLES);
536  glNormal3f(vN->x[0],vN->x[1],vN->x[2]);
537  glVertex3d(v00->x[0],v00->x[1],v00->x[2]);
538  glVertex3d(v11->x[0],v11->x[1],v11->x[2]);
539  glVertex3d(v10->x[0],v10->x[1],v10->x[2]);
540  glEnd();
541  glPopMatrix();//Particle
542 }
543 void ElPoly::DrTria(Vettore *v00,Vettore *v01,Vettore *v11,Vettore *vN){
544  vN->NormalSurf(v00,v01,v11);
545  vN->Normalize();
546  glPushMatrix();//Particle
547  glBegin(GL_TRIANGLES);
548  glNormal3f(vN->x[0],vN->x[1],vN->x[2]);
549  glVertex3d(v00->x[0],v00->x[1],v00->x[2]);
550  glVertex3d(v01->x[0],v01->x[1],v01->x[2]);
551  glVertex3d(v11->x[0],v11->x[1],v11->x[2]);
552  glEnd();
553  glPopMatrix();//Particle
554 }
556  glPushMatrix();//Particle
557  glBegin(GL_LINES);
558  glVertex3d(v00->x[0],v00->x[1],v00->x[2]);
559  glVertex3d(v01->x[0],v01->x[1],v01->x[2]);
560  glEnd();
561  glBegin(GL_LINES);
562  glVertex3d(v01->x[0],v01->x[1],v01->x[2]);
563  glVertex3d(v11->x[0],v11->x[1],v11->x[2]);
564  glEnd();
565  glBegin(GL_LINES);
566  glVertex3d(v11->x[0],v11->x[1],v11->x[2]);
567  glVertex3d(v00->x[0],v00->x[1],v00->x[2]);
568  glEnd();
569  glPopMatrix();//Particle
570 }
571 void ElPoly::DrSmooth(double *Plot,int NSample,double Min,double Max){
572  if(Dr->IfMaterial){
573  glEnable(GL_LIGHTING);
574  glEnable( GL_LIGHT0 );
575  }
576  else
577  glDisable(GL_LIGHTING);
578  double FactC1 = 1./(double)NSample*InvScaleUn*pEdge(CLat1);
579  double FactC2 = 1./(double)NSample*InvScaleUn*pEdge(CLat2);
580  double FactCN = InvScaleUn*ScaleFact;
581  GLfloat Color[4];
582  Vettore v00(3);
583  Vettore v01(3);
584  Vettore v11(3);
585  Vettore v10(3);
586  Vettore vN(3);
587  Min *= InvScaleUn;
588  Max *= InvScaleUn;
589  double Delta = 1./(Max-Min);
590  for(int s=0;s<NSample-1;s++){
591  for(int ss=0;ss<NSample-1;ss++){
592  v00.x[CLat1]=(s+.5)*FactC1;
593  v00.x[CLat2]=(ss+.5)*FactC2;
594  v00.x[CNorm]=Plot[s*NSample+ss]*FactCN;
595  v01.x[CLat1]=(s+1.5)*FactC1;
596  v01.x[CLat2]=(ss+.5)*FactC2;
597  v01.x[CNorm]=Plot[(s+1)*NSample+ss]*FactCN;
598  v11.x[CLat1]=(s+1.5)*FactC1;
599  v11.x[CLat2]=(ss+1.5)*FactC2;
600  v11.x[CNorm]=Plot[(s+1)*NSample+(ss+1)]*FactCN;
601  v10.x[CLat1]=(s+.5)*FactC1;
602  v10.x[CLat2]=(ss+1.5)*FactC2;
603  v10.x[CNorm]=Plot[s*NSample+ss+1]*FactCN;
604  //glColor4f(.4,.4,.6,1.);
605  //DrDoubleTria(&v00,&v01,&v11,&v10,&vN);
606  glBegin(GL_POLYGON);
607  vN.NormalSurf(&v00,&v01,&v11);
608  vN.Normalize();
609  glNormal3f(vN.x[0],vN.x[1],vN.x[2]);
610  double Depth = (v00.x[CNorm]-Min)*Delta*Saturation+ExtParam;
611  Dr->DepthMap(Depth,Color);
612  glMaterialfv(GL_FRONT_AND_BACK,GL_AMBIENT,Color);
613  glMaterialfv(GL_FRONT_AND_BACK,GL_DIFFUSE,Color);
614  glColor4fv(Color);
615  glVertex3d(v00.x[0],v00.x[1],v00.x[2]);
616  Depth = (v01.x[CNorm]-Min)*Delta*Saturation+ExtParam;
617  Dr->DepthMap(Depth,Color);
618  glMaterialfv(GL_FRONT_AND_BACK,GL_AMBIENT,Color);
619  glMaterialfv(GL_FRONT_AND_BACK,GL_DIFFUSE,Color);
620  glColor4fv(Color);
621  glVertex3d(v01.x[0],v01.x[1],v01.x[2]);
622  Depth = (v11.x[CNorm]-Min)*Delta*Saturation+ExtParam;
623  Dr->DepthMap(Depth,Color);
624  glMaterialfv(GL_FRONT_AND_BACK,GL_AMBIENT,Color);
625  glMaterialfv(GL_FRONT_AND_BACK,GL_DIFFUSE,Color);
626  glColor4fv(Color);
627  glVertex3d(v11.x[0],v11.x[1],v11.x[2]);
628  Depth = (v10.x[CNorm]-Min)*Delta*Saturation+ExtParam;
629  Dr->DepthMap(Depth,Color);
630  glMaterialfv(GL_FRONT_AND_BACK,GL_AMBIENT,Color);
631  glMaterialfv(GL_FRONT_AND_BACK,GL_DIFFUSE,Color);
632  glColor4fv(Color);
633  glVertex3d(v10.x[0],v10.x[1],v10.x[2]);
634  glEnd();
635  }
636  }
637  glDisable(GL_LIGHTING);
638 }
639 void ElPoly::DrSample(int NSample){
640  NSample *= 4;
641  MOMENTI m1 = SampleSurfaceMem(NSample);
642  printf("plane min=%lf av=%lf max=%lf\n",m1.Min,m1.Uno,m1.Max);
643  glDeleteLists(Dr->Particles,1);
644  Dr->Particles = glGenLists(1);
645  glNewList(Dr->Particles,GL_COMPILE);
646  DrSmooth(PlotMem,NSample,.9*m1.Min,.9*m1.Max);
647  glEndList();
648 }
650  glDeleteLists(Dr->Particles,1);
651  Dr->Particles = glGenLists(1);
652  glNewList(Dr->Particles,GL_COMPILE);
653  glDisable(GL_NORMALIZE);
654  double Offset = 1./(double)NEdge;
655  Vettore v1(3);
656  Vettore v2(3);
657  Vettore v3(3);
658  Vettore v4(3);
659  Vettore vN(3);
660  glFrontFace(GL_CCW);
661  //glPolygonMode(GL_FRONT_AND_BACK, GL_FILL);
662  //glShadeModel(GL_SMOOTH);
663  for(int p=0,c=0,t=0;p<pNPart();p++){
664  if( Ln[p].NLink != 3 ) continue;
665  if(c == pChain(p)) continue;
666  if( NPType != BEAD_EVERY && pType(p) != NPType && pType(p) != BEAD_NANO) continue;
667  c = pChain(p);
668  double Red = pVel(p,0)*Saturation;// + pPos(p,2)*Sat;
669  double Green = pVel(p,1)*Saturation;
670  double Blue = pVel(p,2)*Saturation;
671  double Alpha = 1.0;//7.+(Green+Blue);
672  int l1 = Ln[p].Link[0];
673  //if(l1 > Gen->NPart) continue;
674  int l2 = Ln[p].Link[2];
675  //if(l2 > Gen->NPart) continue;
676  int l3 = Ln[p].Link[1];
677  // Green = 0.; Blue = 0.;
678  glColor4f(Red,Green,Blue,Alpha);
679  v1.x[0]=pPos(p,0)*InvScaleUn;
680  v1.x[1]=pPos(p,1)*InvScaleUn;
681  v1.x[2]=pPos(p,2)*InvScaleUn*ScaleFact;
682  v2.x[0]=pPos(l1,0)*InvScaleUn;
683  v2.x[1]=pPos(l1,1)*InvScaleUn;
684  v2.x[2]=pPos(l1,2)*InvScaleUn*ScaleFact;
685  v3.x[0]=pPos(l3,0)*InvScaleUn;
686  v3.x[1]=pPos(l3,1)*InvScaleUn;
687  v3.x[2]=pPos(l3,2)*InvScaleUn*ScaleFact;
688  v4.x[0]=pPos(l2,0)*InvScaleUn;
689  v4.x[1]=pPos(l2,1)*InvScaleUn;
690  v4.x[2]=pPos(l2,2)*InvScaleUn*ScaleFact;
691  vN.NormalSurf(&v1,&v2,&v3);
692  vN.Rescale(.2);
693  glPushMatrix();//Particle
694  glBegin(GL_QUADS);
695  glNormal3f(vN.x[0],vN.x[1],-vN.x[2]);
696  glVertex3d(v4.x[0],v4.x[1],v4.x[2]);
697  glVertex3d(v3.x[0],v3.x[1],v3.x[2]);
698  glVertex3d(v2.x[0],v2.x[1],v2.x[2]);
699  glVertex3d(v1.x[0],v1.x[1],v1.x[2]);
700  glEnd();
701  glPopMatrix();//Particle
702  }
703  glEndList();
704 }
706  int NSample = 24;
707  double Threshold = 0.;
708  int NLevel = 4;
709  double **Plot = (double **) calloc(NLevel,sizeof(double));
710  for(int l=0;l<NLevel;l++){
711  Plot[l] = (double *)calloc(NSample*NSample,sizeof(double));
712  }
713  Stalk(NSample,NLevel,Plot,Threshold);
714  //Drawing
715  glColorMaterial(GL_FRONT_AND_BACK,GL_EMISSION);
716  glPolygonMode(GL_FRONT_AND_BACK, GL_FILL);
717  double FactC1 = 1./(double)NSample*InvScaleUn*pEdge(CLat1);
718  double FactC2 = 1./(double)NSample*InvScaleUn*pEdge(CLat2);
719  double FactCN = InvScaleUn*ScaleFact;
720  Vettore v00(3);
721  Vettore v01(3);
722  Vettore v11(3);
723  Vettore v10(3);
724  Vettore vN(3);
725  glDeleteLists(Dr->Particles,1);
726  Dr->Particles = glGenLists(1);
727  glNewList(Dr->Particles,GL_COMPILE);
728  for(int s=0;s<NSample-1;s++){
729  for(int ss=0;ss<NSample-1;ss++){
730  for(int l=0;l<NLevel;l++){
731  if(l==0 || l == 3) continue;
732  glColor4f((double)l*.25,.5,.8,1.);
733  v00.x[CLat1]=(s+.5)*FactC1;
734  v00.x[CLat2]=(ss+.5)*FactC2;
735  v00.x[CNorm]=Plot[l][s*NSample+ss]*FactCN;
736  v01.x[CLat1]=(s+1.5)*FactC1;
737  v01.x[CLat2]=(ss+.5)*FactC2;
738  v01.x[CNorm]=Plot[l][(s+1)*NSample+ss]*FactCN;
739  v11.x[CLat1]=(s+1.5)*FactC1;
740  v11.x[CLat2]=(ss+1.5)*FactC2;
741  v11.x[CNorm]=Plot[l][(s+1)*NSample+ss+1]*FactCN;
742  v10.x[CLat1]=(s+.5)*FactC1;
743  v10.x[CLat2]=(ss+1.5)*FactC2;
744  v10.x[CNorm]=Plot[l][s*NSample+ss+1]*FactCN;
745  if(Plot[l][(s)*NSample+ss] <= 0.){
746  if(l==0)
747  v00.x[CNorm]=Plot[3][(s)*NSample+ss+1]*FactCN;
748  if(l==1)
749  v00.x[CNorm]=Plot[2][(s)*NSample+ss+1]*FactCN;
750  if(l==2)
751  v00.x[CNorm]=Plot[1][(s)*NSample+ss+1]*FactCN;
752  if(l==3)
753  v00.x[CNorm]=Plot[0][(s)*NSample+ss+1]*FactCN;
754  if(v00.x[CNorm] <= 0.) continue;
755  }
756  if(Plot[l][(s+1)*NSample+ss] <= 0.){
757  if(l==0)
758  v01.x[CNorm]=Plot[3][(s)*NSample+ss]*FactCN;
759  if(l==1)
760  v01.x[CNorm]=Plot[2][(s)*NSample+ss]*FactCN;
761  if(l==2)
762  v01.x[CNorm]=Plot[1][(s)*NSample+ss]*FactCN;
763  if(l==3)
764  v01.x[CNorm]=Plot[0][(s)*NSample+ss]*FactCN;
765  if(v01.x[CNorm] <= 0.) continue;
766  }
767  if(Plot[l][s*NSample+ss+1] <= 0.){
768  if(l==0)
769  v10.x[CNorm]=Plot[3][(s)*NSample+ss]*FactCN;
770  if(l==1)
771  v10.x[CNorm]=Plot[2][(s)*NSample+ss]*FactCN;
772  if(l==2)
773  v10.x[CNorm]=Plot[1][(s)*NSample+ss]*FactCN;
774  if(l==3)
775  v10.x[CNorm]=Plot[0][(s)*NSample+ss]*FactCN;
776  if(v10.x[CNorm] <= 0.) continue;
777  }
778  if(Plot[l][(s+1)*NSample+(ss+1)] <= 0.){
779  if(l==0)
780  v11.x[CNorm]=Plot[3][(s)*NSample+ss+1]*FactCN;
781  if(l==1)
782  v11.x[CNorm]=Plot[2][(s)*NSample+ss+1]*FactCN;
783  if(l==2)
784  v11.x[CNorm]=Plot[1][(s)*NSample+ss+1]*FactCN;
785  if(l==3)
786  v11.x[CNorm]=Plot[0][(s)*NSample+ss+1]*FactCN;
787  if(v11.x[CNorm] <= 0.) continue;
788  }
789  DrDoubleTria(&v00,&v01,&v11,&v10,&vN);
790  }
791  }
792  }
793  glEndList();
794  //Freeing
795  for(int l=0;l<NLevel;l++){
796  free(Plot[l]);
797  }
798  free(Plot);
799 }
801  int NSample = 32;
802  double Unity = 1./(double)NSample*pEdge(0)*InvScaleUn;
803  double *Plot = (double *) calloc(SQR(NSample),sizeof(double));
804  double *PlotR = (double *) calloc(SQR(NSample),sizeof(double));
805  SampleSurface(Plot,NSample,NChType);
806  //InterBSpline2D(Plot,PlotR,NSample,NSample);
807  Spettro2d(Plot,NSample);
808  double Min = 0.;
809  double Max = 1000000.;
810  for(int s=0;s<SQR(NSample);s++){
811  Plot[s] += .5*pEdge(CNorm)/ScaleFact;
812  if(Max < Plot[s]) Max = Plot[s];
813  if(Min > Plot[s]) Min = Plot[s];
814  }
815  //Drawing
816  glDeleteLists(Dr->Particles,1);
817  Dr->Particles = glGenLists(1);
818  glNewList(Dr->Particles,GL_COMPILE);
819  DrSmooth(Plot,NSample,Min,Max);
820  glEndList();
821  //Freeing
822  free(PlotR);
823  free(Plot);
824  return;
825  double *Points = (double *)calloc(NSample,sizeof(double));
826  Spettro2d(Points,NSample,NChType);
827  double Max1=0.,Max2=0.;
828  for(int s=0;s<NSample;s++){
829  if(Max1 < Points[s]){
830  Max2 = Max1;
831  Max1 = Points[s];
832  continue;
833  }
834  if(Max2 < Points[s]){
835  Max2 = Points[s];
836  continue;
837  }
838  }
839  for(int s=0;s<NSample;s++){
840  Points[s] /= Max2;
841  }
842  //SampleSurface(Plot,NSample,2,CHAIN_DOWN);
843  for(int s=NSample/2;s<NSample-1;s++){
844  glPushMatrix();
845  glBegin(GL_LINES);
846  glColor4f(1.,0.,1.,1.);
847  glVertex3d((2*s-NSample)*Unity,0.,Points[s]);
848  glVertex3d((2*(s+1)-NSample)*Unity,0.,Points[s+1]);
849  glEnd();
850  glPopMatrix();
851  }
852  glEndList();
853  free(Points);
854  free(Plot);
855 }
857  int NSample = 20;
858  double Unity= 1./(double)NSample*pEdge(0)*InvScaleUn;
859  SPLINE Weight;
860  Weight.a0 = 1.; Weight.a1 = 0.; Weight.a2=1.;Weight.a3=0.;Weight.a4=0.;
861  Matrice *Surface= new Matrice(NSample,NSample);
862  SampleSurface(Surface,NSample,NChType);
863  Matrice *Resp = new Matrice(NSample,NSample);
864  SpatialDerivative(Surface,Resp,Weight,NSample);
865  glDeleteLists(Dr->Particles,1);
866  Dr->Particles = glGenLists(1);
867  glNewList(Dr->Particles,GL_COMPILE);
868  for(int s=1;s<NSample-2;s++){
869  for(int ss=1;ss<NSample-2;ss++){
870  //printf("Plot[%d][%d]= %lf\n",s,ss,Plot[s][ss]);
871  glPushMatrix();//Particle
872  glColor4f(0.,1.,1.,1.);
873  glBegin(GL_TRIANGLES);
874  glNormal3f(0.,0.,1.);
875  glVertex3d((s+.5)*Unity,(ss+.5)*Unity,Resp->Val(s,ss)*InvScaleUn);
876  glVertex3d((s+1.5)*Unity,(ss+.5)*Unity,Resp->Val(s+1,ss)*InvScaleUn);
877  glVertex3d((s+1.5)*Unity,(ss+1.5)*Unity,Resp->Val(s+1,ss+1)*InvScaleUn);
878  glEnd();
879  glPopMatrix();//Particle
880  glPushMatrix();//Particle
881  glColor4f(0.5,.0,0.5,1.);
882  glBegin(GL_TRIANGLES);
883  glNormal3f(0.,0.,1.);
884  glVertex3d((s+.5)*Unity,(ss+.5)*Unity,Resp->Val(s,ss)*InvScaleUn);
885  glVertex3d((s+1.5)*Unity,(ss+1.5)*Unity,Resp->Val(s+1,ss+1)*InvScaleUn);
886  glVertex3d((s+.5)*Unity,(ss+1.5)*Unity,Resp->Val(s,ss+1)*InvScaleUn);
887  glEnd();
888  glPopMatrix();//Particle
889  }
890  }
891  delete Surface;
892  delete Resp;
893  glEndList();
894 }
895 void ElPoly::DrIsoipse(int Values,int NIsoipse,int CoordN){
896  if(CoordN < 0 || CoordN > 3) return;
897  glDeleteLists(Dr->Particles,1);
898  Dr->Particles = glGenLists(1);
899  glNewList(Dr->Particles,GL_COMPILE);
900  double InvIsoipse = 1./(double)NIsoipse;
901  int Coord1 = (CoordN+1)%3;
902  int Coord2 = (CoordN+2)%3;
903  int *Plot = (int *) calloc(Values*Values*NIsoipse,sizeof(int));
904  for(int p=0;p<pNPart();p++){
905  int v = (int) (pPos(p,Coord1)*NIsoipse/pEdge(Coord1));
906  if( v < 0 || v > Values) continue;
907  int vv = (int) (pPos(p,Coord2)*NIsoipse/pEdge(Coord2));
908  if( vv < 0 || vv > Values) continue;
909  int i = (int) (pPos(p,CoordN)*NIsoipse/pEdge(CoordN));
910  if( i < 0 || i > NIsoipse) continue;
911  Plot[(v*Values + vv)*NIsoipse + i] = p;
912  }
913  for(int i=0;i<NIsoipse;i++){
914  glBegin(GL_POLYGON);
915  glNormal3f(0.,0.,1.);
916  double Zed = i*pEdge(CoordN)*InvIsoipse*InvScaleUn*ScaleFact;
917  for(int p=0;p<pNPart();p++){
918  int v = (int) (pPos(p,CoordN)*NIsoipse/pEdge(CoordN));
919  if( v != i) continue;
920  for(int v=0;v<Values-1;v++){
921  int pMin = 0;
922  int pMax = 0;
923  for(int vv=0;vv<Values;vv++){
924  if(Plot[(v*Values+vv)*NIsoipse+i] == 0 && Plot[(v*Values+vv+1)*NIsoipse+i] != 0)
925  pMin = Plot[(v*Values+vv+1)*NIsoipse+i];
926  if(Plot[(v*Values+vv)*NIsoipse+i] != 0 && Plot[(v*Values+vv+1)*NIsoipse+i] == 0){
927  pMax = Plot[(v*Values+vv)*NIsoipse+i];
928  break;
929  }
930  }
931  glVertex3d(pPos(pMin,Coord1)*InvScaleUn,pPos(pMin,Coord2)*InvScaleUn,Zed);
932  glVertex3d(pPos(pMax,Coord1)*InvScaleUn,pPos(pMax,Coord2)*InvScaleUn,Zed);
933  }
934  glEnd();
935  }
936  }
937  free(Plot);
938  glEndList();
939 }
941  int NSample = 40;
942  PART *Pn = (PART *)calloc(NSample*NSample,sizeof(PART));
943 }
945  int NDim = 2;
946  int NCell = pNChain();
947  int Batch = 1000;
948  int DontCreate = 4;
949  int SampleType = 1;//Halton
950  int NSample = 1000;
951  int NItMax = 4;
952  int NItFix = 1;
953  int Seme = 4532643;
954  int NIt = 0;
955  double ItDiff = 0.;
956  double Energy = 0.;
957  double *Pos = (double *)calloc(NCell*NDim,sizeof(double));
958  for(int c=0;c<NCell;c++){
959  for(int d=0;d<NDim;d++){
960  Pos[c*NDim+d] = Ch[c].Pos[d];
961  }
962  }
963  //cvt(NDim,NCell,Batch,DontCreate,SampleType,NSample,NItMax,NItFix,&Seme,Pos,&NIt,&ItDiff,&Energy);
964  Hexagon = Dr->DefHexagon();
965  glDeleteLists(Dr->Particles,1);
966  Dr->Particles = glGenLists(1);
967  glNewList(Dr->Particles,GL_COMPILE);
968  for(int c=0;c<NCell;c++){
969  glPushMatrix();
970  glColor4f(.2,.5,.3,1.);
971  glBegin(GL_POINTS);
972  glNormal3f(0.,0.,1.);
973  glVertex3d(Pos[c*NDim+0]*InvScaleUn,
974  Pos[c*NDim+1]*InvScaleUn,.5*pEdge(2)*InvScaleUn);
975  glEnd();
976  //glCallList(Hexagon);
977  glPopMatrix();
978  glColor4f(.8,.2,.5,1.);
979  glPushMatrix();
980  glTranslated(Ch[c].Pos[0]*InvScaleUn,
981  Ch[c].Pos[1]*InvScaleUn,.5*pEdge(2)*InvScaleUn);
982  glCallList(Hexagon);
983  glPopMatrix();
984  }
985  glEndList();
986 }
988  DrVoronoi();
989  int Vertex = 6;
990  int **Triangle = (int **)calloc(Vertex+1,sizeof(int));
991  for(int v=0;v<Vertex+1;v++){
992  *(Triangle+v) = (int *)calloc(pNChain(),sizeof(int));
993  }
994  //Arrange(Triangle,Vertex);
995  glDeleteLists(Dr->Particles,1);
996  Dr->Particles = glGenLists(1);
997  glNewList(Dr->Particles,GL_COMPILE);
998  for(int c=0;c<pNChain();c++){
999  int Chc = Ch[c].Type;
1000  if(!CHAIN_IF_TYPE(Chc,NChType)) continue;
1001  // glBegin(GL_POLYGON);
1002  int l = Triangle[Vertex-1][c];
1003  int l1 = Triangle[0][c];
1004  glPushMatrix();//Particle
1005  // double PosX1 = .333*(Ch[l].Pos[0]+Ch[c].Pos[0]+Ch[l1].Pos[0])*InvScaleUn;
1006  // double PosY1 = .333*(Ch[l].Pos[1]+Ch[c].Pos[1]+Ch[l1].Pos[1])*InvScaleUn;
1007  // double PosZ1 = .333*(Ch[l].Pos[2]+Ch[c].Pos[2]+Ch[l1].Pos[2])*InvScaleUn;
1008  for(int v=0;v<Vertex;v++){
1009  l = Triangle[v][c];
1010  if(v<Vertex-1) l1 = Triangle[v+1][c];
1011  else l1 = Triangle[0][c];
1012  // double PosX = .333*(Ch[l].Pos[0]+Ch[c].Pos[0]+Ch[l1].Pos[0])*InvScaleUn;
1013  // double PosY = .333*(Ch[l].Pos[1]+Ch[c].Pos[1]+Ch[l1].Pos[1])*InvScaleUn;
1014  // double PosZ = .333*(Ch[l].Pos[2]+Ch[c].Pos[2]+Ch[l1].Pos[2])*InvScaleUn;
1015  double Shift[3] = {0.,0.,0.};
1016  for(int d=0;d<3;d++) Shift[d] = -floor((Ch[c].Pos[d]-Ch[l].Pos[d])/pEdge(d)+.5)*pEdge(d);
1017  glPushMatrix();//Line
1018  glBegin(GL_LINES);
1019  glNormal3f(0.,0.,1.);
1020  glVertex3f(Ch[c].Pos[0]*InvScaleUn,Ch[c].Pos[1]*InvScaleUn,Ch[c].Pos[2]*InvScaleUn*ScaleFact);
1021  glVertex3f((Ch[l].Pos[0]-Shift[0])*InvScaleUn,(Ch[l].Pos[1]-Shift[1])*InvScaleUn,(Ch[l].Pos[2]-Shift[2])*InvScaleUn*ScaleFact);
1022  glEnd();
1023  //PosX1 = PosX;PosY1=PosY;PosZ1=PosZ;
1024  glPopMatrix();//Line
1025  }
1026  // printf("\n");
1027  // printf("(%lf %lf %lf) \n",Ch[c].Pos[0],Ch[c].Pos[1],Ch[c].Pos[2]);
1028  glColor4f(1.,1.,0.,1.);
1029  //glEnd();
1030  glPopMatrix();//Particle
1031  glPushMatrix();
1032  glTranslated(Ch[c].Pos[0]*InvScaleUn,
1033  Ch[c].Pos[1]*InvScaleUn,
1034  Ch[c].Pos[2]*InvScaleUn*ScaleFact);
1035  glColor4f(1.,0.,0.,1.);
1036  glCallList(Point);
1037  // glutSolidSphere(.05,20,20);
1038  glPopMatrix();
1039  glPushMatrix();//Info
1040  glColor3f(1.0,1.0,1.0);
1041  // glTranslated(Ch[c].Pos[0]*InvScaleUn,Ch[c].Pos[1]*InvScaleUn,Ch[c].Pos[2]*InvScaleUn);
1042  glRasterPos3f(Ch[c].Pos[0]*InvScaleUn,
1043  Ch[c].Pos[1]*InvScaleUn,
1044  Ch[c].Pos[2]*InvScaleUn*ScaleFact);
1045  char *Numero = (char*)calloc(20,sizeof(char));
1046  sprintf(Numero,"%d",c);//CalcnPos(Ch[c].Pos));
1047  for(int i=0;i<strlen(Numero);i++)
1048  glutBitmapCharacter(GLUT_BITMAP_TIMES_ROMAN_24,Numero[i]);
1049  free(Numero);
1050  glPopMatrix();//Info
1051  }
1052  glEndList();
1053  free(Triangle);
1054  // GLuint *Surface = (GLuint *)malloc(3*Gen->NChain*sizeof(GLuint));
1055  // for(int c=0,cc=0;c<Gen->NChain;c++){
1056  // for(int d=0;d<3;d++){
1057  // Surface[3*c+d] = (GLuint) Ch[c].Pos[d];
1058  // }
1059  // printf("(%u %u %u) (%lf %lf %lf) \n",Surface[3*c+0],Surface[3*c+1],Surface[3*c+2],Ch[c].Pos[0],Ch[c].Pos[1],Ch[c].Pos[2]);
1060  // }
1061  //glEnable(GL_VERTEX_ARRAY);
1062  // glDrawElements(GL_TRIANGLES,Gen->NChain,GL_UNSIGNED_INT,Surface);
1063  //free(Surface);
1064 }
1066  double Dx = sqrt(pEdge(CLat1)*pEdge(CLat2)/(double)pNPart());
1067  int NLat1 = (int)(pEdge(CLat1)/Dx);
1068  int NLat2 = (int)(pEdge(CLat2)/Dx);
1069  printf("Dx %lf NLat1 %d NLat2 %d %d\n",Dx,NLat1,NLat2,SQR(NLat1) + SQR(NLat2));
1070  if( NLat1*NLat2 != pNPart()) return;
1071  glDeleteLists(Dr->Particles,1);
1072  Dr->Particles = glGenLists(1);
1073  glNewList(Dr->Particles,GL_COMPILE);
1074  Vettore v1(3);
1075  Vettore v2(3);
1076  Vettore v3(3);
1077  Vettore v4(3);
1078  Vettore vN(3);
1079  printf("Ciccia\n");
1080  for(int nx = 0;nx<NLat1-2;nx+=2){
1081  for(int ny = 0;ny<NLat2-2;ny+=2){
1082  int l1 = nx*NLat2 + ny;
1083  int l2 = nx*NLat2 + ny + 1;
1084  int l3 = (nx+1)*NLat2 + ny + 1;
1085  int l4 = (nx+1)*NLat2 + ny;
1086  v1.x[0]=pPos(l1,0)*InvScaleUn;
1087  v1.x[1]=pPos(l1,1)*InvScaleUn;
1088  v1.x[2]=pPos(l1,2)*InvScaleUn*ScaleFact;
1089  v2.x[0]=pPos(l2,0)*InvScaleUn;
1090  v2.x[1]=pPos(l2,1)*InvScaleUn;
1091  v2.x[2]=pPos(l2,2)*InvScaleUn*ScaleFact;
1092  v3.x[0]=pPos(l3,0)*InvScaleUn;
1093  v3.x[1]=pPos(l3,1)*InvScaleUn;
1094  v3.x[2]=pPos(l3,2)*InvScaleUn*ScaleFact;
1095  v4.x[0]=pPos(l4,0)*InvScaleUn;
1096  v4.x[1]=pPos(l4,1)*InvScaleUn;
1097  v4.x[2]=pPos(l4,2)*InvScaleUn*ScaleFact;
1098  vN.NormalSurf(&v1,&v2,&v3);
1099  DrDoubleTria(&v1,&v2,&v3,&v4,&vN);
1100  }
1101  }
1102  glEndList();
1103 }
1104 void ElPoly::DrIsolevel(int NSample,double IsoLevel){
1105  NSample = 20;
1106  double VolEl = pVol()/(double)CUB(NSample);
1107  double *Plot = (double *)calloc(CUBE(NSample),sizeof(double));
1108  double Min = 0.;
1109  double Max = 0.;
1110  VAR_TRIANGLE *Tri = NULL;
1111  for(int p=0;p<pNPart();p++){
1112  //if(Pm[p].Typ != 1) continue;
1113  double Posx = pPos(p,0) - floor(pPos(p,0)*pInvEdge(0))*pEdge(0);
1114  double Posy = pPos(p,1) - floor(pPos(p,1)*pInvEdge(1))*pEdge(1);
1115  double Posz = pPos(p,2) - floor(pPos(p,2)*pInvEdge(2))*pEdge(2);
1116  int sx = (int)(Posx*pInvEdge(0)*NSample);
1117  int sy = (int)(Posy*pInvEdge(1)*NSample);
1118  int sz = (int)(Posz*pInvEdge(2)*NSample);
1119  int sTot = (sx*NSample+sy)*NSample+sz;
1120  Plot[sTot] += VolEl;
1121  if(Max < Plot[sTot]) Max = Plot[sTot];
1122  if(Min > Plot[sTot]) Min = Plot[sTot];
1123  }
1124  //IsoLevel = .1*Max;
1125  printf("Min %lf Max %lf Isolevel %lf\n",Min,Max,IsoLevel);
1126  int NTri = 0;
1127  Tri = MarchingCubes(Plot,NSample,IsoLevel,&NTri);
1128  free(Plot);
1129  int NVertex = CUBE(2*NSample-1);
1130  double *WeightL = (double *) calloc(NVertex,sizeof(double));
1131  double MaxW = 1.;
1132  for(int t=0;t<NTri;t++){
1133  int vRef = Tri[t].v[0];
1134  //printf("%d %lf\n",t,WeightL[vRef]);
1135  }
1136  glDeleteLists(Dr->Particles,1);
1137  Dr->Particles = glGenLists(1);
1138  glNewList(Dr->Particles,GL_COMPILE);
1139  for(int n=0;n<pNNano();n++) DrawNano(n);
1140  for(int t=0;t<NTri;t++){
1141  //if(Tri[t].p[0].x[2]> Nano->Pos[2] + 3.||Tri[t].p[0].x[2] < Nano->Pos[2] - 3.) continue;
1142  int vRef = Tri[t].v[0];
1143  double Green = WeightL[vRef]*MaxW + .3 + .4*Mat->Casuale();
1144  glColor4f(0.1,Green,0.2,1.);
1145  //glColor4f(HueUp[p].r,HueUp[p].g,HueUp[p].b,HueUp[p].a);
1146  glBegin(GL_TRIANGLES);
1147  glNormal3f(Tri[t].n[0].x[0]*InvScaleUn,
1148  Tri[t].n[0].x[1]*InvScaleUn,
1149  Tri[t].n[0].x[2]*InvScaleUn*ScaleFact);
1150  glVertex3f(Tri[t].p[0].x[0]*InvScaleUn,
1151  Tri[t].p[0].x[1]*InvScaleUn,
1152  Tri[t].p[0].x[2]*InvScaleUn*ScaleFact);
1153  glNormal3f(Tri[t].n[1].x[0]*InvScaleUn,
1154  Tri[t].n[1].x[1]*InvScaleUn,
1155  Tri[t].n[1].x[2]*InvScaleUn*ScaleFact);
1156  glVertex3f(Tri[t].p[1].x[0]*InvScaleUn,
1157  Tri[t].p[1].x[1]*InvScaleUn,
1158  Tri[t].p[1].x[2]*InvScaleUn*ScaleFact);
1159  glNormal3f(Tri[t].n[2].x[0]*InvScaleUn,
1160  Tri[t].n[2].x[1]*InvScaleUn,
1161  Tri[t].n[2].x[2]*InvScaleUn*ScaleFact);
1162  glVertex3f(Tri[t].p[2].x[0]*InvScaleUn,
1163  Tri[t].p[2].x[1]*InvScaleUn,
1164  Tri[t].p[2].x[2]*InvScaleUn*ScaleFact);
1165  glEnd();
1166  if(IfLine){
1167  glColor4f(0.1,.3,1.,1.);
1168  glBegin(GL_LINES);
1169  glVertex3d(Tri[t].c.x[0]*InvScaleUn,
1170  Tri[t].c.x[1]*InvScaleUn,
1171  Tri[t].c.x[2]*InvScaleUn*ScaleFact);
1172  glVertex3d((Tri[t].c.x[0]+Tri[t].n[0].x[0])*InvScaleUn,
1173  (Tri[t].c.x[1]+Tri[t].n[0].x[1])*InvScaleUn,
1174  (Tri[t].c.x[2]+Tri[t].n[0].x[2])*InvScaleUn*ScaleFact);
1175  glEnd();
1176  }
1177  }
1178  glEndList();
1179  free(WeightL);
1180 }
1181 #endif //USE_GL
void DrDensity()
(rad,norm,dens) visualisation
XYZ n[3]
Normal to the vertices.
Definition: VarData.h:506
double a4
a0 + a1*x + a2*x^2 + a3*x^3 + a4^4
CHAIN * Ch
Information on all chains.
Definition: VarData.h:1050
void DrSmooth(double *Plot, int NSample, double Min, double Max)
Draw the surface in triangles.
void DrStalk()
Visualize the surface from Stalk.xvl.
double InvScaleUn
Define the shrink factor between the box edges.
Definition: ElPoly.h:457
double x[3]
Cartesian coordinates.
Definition: VarData.h:486
double NanoDist2(double *Pos, int n)
Pointer to a generic function.
Definition: VarData.h:764
int Type
Type of the chain (see list CHAIN_)
Definition: VarData.h:246
double pVel(int p, int d)
Return the velocity.
LINKS * Ln
Array of linking between the particles.
Definition: VarData.h:1048
Draw provides the basic configuration of the openGL libraries used in every derived program...
Definition: Draw.h:15
void DrIsoipse(int NBin, int NIsoipse, int CoordN)
Draws all the same quotes surfaces.
void DrQuad()
Call a square in every particle position.
void DrCreateStalk()
Visualize the surface calculated in Stalk()
double Normalize()
Normlizes a Vettore.
void DrIsolevel(int NSample, double IsoLevel)
Defines the triangles at the boundaries of the density close to the IsoLevel value.
void DrPolygon()
Every particle is traeted as a vertex.
int pChain(int p)
Return the chain.
Geometrical operations on vectors.
Definition: MatematicaVect.h:9
void DrPotential()
Potential function.
void DrDerivative()
Samples the surface, applies a 2d matrix derivative and draws the surface in triangles.
double Min
Minimum value.
double pInvEdge(int d)
Inverted xyzr edges of the simulation box.
Definition: VarData.h:920
void DrTriaContour(Vettore *v00, Vettore *v01, Vettore *v11)
Draw the contour of three vectors.
int NChType
Type of chain selected.
Definition: VarData.h:1080
void DrNormalPoint(int p, int NEdge)
Draw the normal to a point.
Define a triangle.
Definition: VarData.h:498
MOMENTI SampleSurfaceMem(int NSample)
Allocate and fill PlotMem with the particle average position.
Definition: VarDataEl.cpp:49
int pType(int p)
Return the type.
int CLat2
lateral coordinate
Definition: VarData.h:1078
VAR_TRIANGLE * MarchingCubes(double *Plot, int NSample, double IsoLevel, int *NTri)
Defines the triangles close to the IsoLevel of the 3d density Plot.
void DepthMap(double Val, GLfloat *Color)
Pointer to a generic function.
Definition: Draw.h:113
void DrSample(int NSample)
Samples the surface and draws it in triangles.
double * PlotMem
Particle position/density on the square lattice.
Definition: VarData.h:1056
double * x
Where the data are stored.
double ScaleFact
Normal scaling factor (z zooming)
Definition: ElPoly.h:461
void DrDoubleTria(Vettore *v00, Vettore *v01, Vettore *v11, Vettore *v10, Vettore *vN)
Double triangle defined by four vectors.
void Rescale(double NewLength)
Rescale the total length of the vector.
double pVol()
xyzr edges of the simulation box
Definition: VarData.h:922
int v[3]
Reference for the vectors.
Definition: VarData.h:504
void DrQuad1()
Intensity scheme.
double pEdge(int d)
xyzr edges of the simulation box
Definition: VarData.h:918
double Saturation
Saturation of the color (increase intensity)
Definition: ElPoly.h:388
double a3
a0 + a1*x + a2*x^2 + a3*x^3 + a4^4
Moments of a distribution.
int NPType
Type of particle selected.
Definition: VarData.h:1082
int NEdge
Number of particles per edge.
Definition: VarData.h:1084
double Casuale()
Random uniform number.
double Val(int row)
Returns a value in 1d.
void DrawNano(int n)
Pointer to a generic function.
Definition: ElPoly.h:308
int IfMaterial
Activate the illumination for a specific material.
Definition: Draw.h:240
void DrVoronoi()
Voronoi tassellation.
void DrSurface()
Draws the surface like a sheet using the chains position.
double ExtParam
External parameter (e.g. for MarchingCubes)
Definition: ElPoly.h:459
int pNNano()
Number of nanoparticles.
double a2
a0 + a1*x + a2*x^2 + a3*x^3 + a4^4
double a1
a0 + a1*x + a2*x^2 + a3*x^3 + a4^4
GLuint Particles
Refers to the list of the total position of the particles which will be generated in another program...
Definition: Draw.h:200
void Stalk(int NSample, int NLevel, double **Plot, double Threshold)
Following the contour of a stalk.
void DrShell()
Draw the outer shell.
int DefHexagon()
Hexagon.
void DrSpectrum()
Calculate the spectrum of the surface and draws it in triangles.
double Pos[4]
xyzr Postion of the chain
Definition: VarData.h:238
void DrInterpSurface()
Boh.
Matrice computes the algebric operations on matrices.
void Spettro2d(double *Points, int NSample, int Type)
1-d spectrum of a surface
Definition: VarDataExp.cpp:216
void Tile()
Cover a regular square grid of points with tiles.
int IfLine
Boh.
Definition: ElPoly.h:435
double Max
Maximum value.
double Uno
First moment.
Matematica * Mat
Implementation of all usefull algorythms.
Definition: VarData.h:527
double pVelMax(int d)
Maximum velocity.
Definition: VarData.h:926
double pPos(int p, int d)
Return back folded position.
int pNChain()
Number of chain.
int CLat1
lateral coordinate
Definition: VarData.h:1076
Information of every particle.
Definition: VarData.h:214
int SpatialDerivative(Matrice *Surface, Matrice *Resp, SPLINE Weight, int NSample)
Perform a spatial derivative on a surface.
Definition: VarDataEl.cpp:222
int Surface(int NBin, int Coord)
Area of a surfuce projected on the coordinate Coord surf.
Definition: ElPolyRepr.cpp:111
double a0
a0 + a1*x + a2*x^2 + a3*x^3 + a4^4
Coefficient of a spline.
void Mult(double Fact)
multiply by a scalar
void DrSquareMesh()
Draw the polygons for a square mesh.
PART * Pm
Particle information of all particle.
Definition: VarData.h:1046
void DrTria(Vettore *v00, Vettore *v01, Vettore *v11, Vettore *vN)
single triangle defined by three vectors
void SampleSurface(double *Plot, int NSample, int Type)
Define a normal coordinate for every patch.
Definition: VarDataEl.cpp:65
XYZ c
Centroid.
Definition: VarData.h:502
void NormalSurf(const Vettore *u, const Vettore *v, const Vettore *w)
Computes the normal to the plane described by u,v,w.
void DrField(int NGrid, double IsoLevel, int nNano, FILE *FWrite)
Draw a scalar field.
int CNorm
Normal coordinate.
Definition: VarData.h:1074
int pNPart()
Number of particle.