implemented bmp writing (doesntwork by now)
[physik/nlsop.git] / nlsop.c
1 /*
2  * nlsop.c 
3  *
4  * this program tries helping to understand the amorphous depuration
5  * and recrystallization of SiCx while ion implantation at temperatures
6  * below 400 degree celsius.
7  * hopefully the program will simulate the stabilization of the
8  * selforganizing lamella structure in the observed behaviour.
9  *
10  * refs: 
11  *  - J. K. N. Lindner. Habilationsschrift, Universitaet Augsburg.
12  *  - Maik Haeberlen. Diplomarbeit, Universitaet Augsburg.
13  */
14
15 #define _GNU_SOURCE
16 #include <stdio.h>
17 #include <stdlib.h>
18 #include <string.h>
19 #include <sys/types.h>
20 #include <sys/stat.h>
21 #include <fcntl.h>
22 #include <unistd.h>
23
24 #include "nlsop.h"
25
26 #include "dfbapi.h"
27 #include "random.h"
28
29 #define MAKE_AMORPH(N) *(N)|=AMORPH
30 #define MAKE_CRYST(N) *(N)&=~AMORPH
31
32 /* test globals - get removed or included in my_info struct later */
33 int amorph_count;
34 int cryst_count;
35
36 int usage(void)
37 {
38  puts("usage:");
39  puts("-h \t\t help");
40  puts("-n \t\t no user interaction");
41  puts("-Z \t\t cryst -> amorph c diffusion in z direction");
42  printf("-a <value> \t slope of nuclear energy loss (default %f)\n",A_EL);
43  printf("-b <value> \t nuclear energy loss offset (default %f)\n",B_EL);
44  printf("-x <value> \t # x cells (default %d)\n",X);
45  printf("-y <value> \t # y cells (default %d)\n",Y);
46  printf("-z <value> \t # z cells (default %d)\n",Z);
47  printf("-s <value> \t steps (default %d)\n",STEPS);
48  printf("-d <value> \t refresh display (default %d)\n",REFRESH);
49  printf("-r <value> \t amorphous influence range (default %d)\n",RANGE);
50  printf("-f <value> \t pressure = <value> * 1/distance^2 (default %f)\n",A_AP);
51  printf("-p <value> \t pressure offset (default %f)\n",B_AP);
52  printf("-F <value> \t proportionality constant between c conc and ability to get amorphous (default %f)\n",A_CP);
53  printf("-A <value> \t slope of linear c distribution (default %f)\n",A_CD);
54  printf("-B <value> \t linear c distribution offset (default %f)\n",B_CD);
55  printf("-D <value> \t diffusion rate from cryst to amorph cells (default %f)\n",DR_AC);
56  printf("-c <value> \t diffusion rate in cryst cells (default %f)\n",DR_CC);
57  printf("-e <value> \t do diffusion every <value> steps (default %d)\n",DIFF_RATE);
58  printf("-W <value> \t write every <value> steps to save file (default %d)\n",RESAVE);
59  puts("-C <file> \t convert file to gnuplot format");
60  puts("-L <file> \t load from file");
61  puts("-S <file> \t save to file");
62  puts("-R <file> \t read from random file");
63  puts("-P <file> \t specify implantatin profile file");
64  
65  return 1;
66 }
67
68 int process_cell(d3_lattice *d3_l,u32 x,u32 y,u32 z,info *my_info)
69 {
70  unsigned char *thiz;
71  int *conc;
72  int i,j;
73  int off;
74  double p;
75
76  thiz=d3_l->status+x+y*d3_l->max_x+z*d3_l->max_x*d3_l->max_y;
77  conc=d3_l->extra+x+y*d3_l->max_x+z*d3_l->max_x*d3_l->max_y;
78  p=my_info->b_ap*URAND_MAX;
79  for(i=-(my_info->range);i<=my_info->range;i++)
80  {
81   for(j=-(my_info->range);j<=my_info->range;j++)
82   {
83    if(!(i==0 && j==0))
84    {
85     off=((x+d3_l->max_x+i)%d3_l->max_x)+((y+d3_l->max_y+j)%d3_l->max_x)*d3_l->max_x+z*d3_l->max_x*d3_l->max_y;
86     if(*(d3_l->status+off)&AMORPH) p+=my_info->a_ap*(*(d3_l->extra+off))*URAND_MAX/(i*i+j*j);
87    } 
88   }
89  }
90  p+=*conc*my_info->a_cp*URAND_MAX;
91  if(!(*thiz&AMORPH))
92  {
93   if(get_rand(URAND_MAX)<=p)
94   {
95    MAKE_AMORPH(thiz);
96    amorph_count++;
97   }
98  } else
99  {
100   /* assume 1-p probability */
101   if(get_rand(URAND_MAX)>p)
102   {
103    MAKE_CRYST(thiz);
104    cryst_count++;
105   }
106  }
107  
108  return 1;
109 }
110
111 int distrib_c(d3_lattice *d3_l,info *my_info,int step,double c_ratio)
112 {
113  u32 x,y,z;
114  int i,j,k,c;
115  int offset,off;
116  int carry;
117
118  /* put one c ion somewhere in the lattice */
119  if(my_info->cc<c_ratio*step)
120  {
121   x=get_rand(d3_l->max_x);
122   y=get_rand(d3_l->max_y);
123   z=get_rand_lgp(d3_l->max_z,my_info->a_cd,my_info->b_cd);
124   *(d3_l->extra+x+y*d3_l->max_x+z*d3_l->max_x*d3_l->max_y)+=1;
125   (my_info->cc)++;
126  }
127
128  if(step%my_info->diff_rate==0)
129  {
130
131  for(i=0;i<d3_l->max_x;i++)
132  {
133   for(j=0;j<d3_l->max_y;j++)
134   {
135    for(k=0;k<d3_l->max_z;k++)
136    {
137     offset=i+j*d3_l->max_x+k*d3_l->max_x*d3_l->max_y;
138     /* case amorph - amorph <-> cryst diffusion */
139     if(*(d3_l->status+offset)&AMORPH)
140     {
141      for(c=-1;c<=1;c++)
142      {
143       if(c!=0)
144       {
145        off=((i+d3_l->max_x+c)%d3_l->max_x)+j*d3_l->max_x+k*d3_l->max_x*d3_l->max_y;
146        carry=0;
147        if(!(*(d3_l->status+off)&AMORPH)) carry=(int)(my_info->dr_ac*(*(d3_l->extra+off)));
148        if(carry!=0)
149        {
150         *(d3_l->extra+offset)+=carry;
151         *(d3_l->extra+off)-=carry;
152        }
153       }
154      }
155      for(c=-1;c<=1;c++)
156      {
157       if(c!=0)
158       {
159        off=i+((j+c+d3_l->max_y)%d3_l->max_y)*d3_l->max_x+k*d3_l->max_x*d3_l->max_y;
160        carry=0;
161        if(!(*(d3_l->status+off)&AMORPH)) carry=(int)(my_info->dr_ac*(*(d3_l->extra+off)));                             
162        if(carry!=0)
163        {
164         *(d3_l->extra+offset)+=carry; 
165         *(d3_l->extra+off)-=carry; 
166        }
167       }
168      }
169      if(my_info->z_diff)
170      {
171       for(c=-1;c<=1;c++)
172       {
173        if(c!=0)
174        {
175         off=i+j*d3_l->max_x+((k+c+d3_l->max_z)%d3_l->max_z)*d3_l->max_x*d3_l->max_y;
176         carry=0;
177         if(!*(d3_l->status+off)&AMORPH) carry=(int)(my_info->dr_ac*(*(d3_l->extra+off)));
178         if(carry!=0)
179         {
180          *(d3_l->extra+off)-=carry;
181          *(d3_l->extra+offset)+=carry;
182         }
183        }
184       }
185      }  
186     } else
187     /* case not amorph - cryst <-> cryst diffusion */
188     {
189      for(c=-1;c<=1;c++) 
190      {
191       if(c!=0)
192       {
193        off=i+((j+c+d3_l->max_y)%d3_l->max_y)*d3_l->max_x+k*d3_l->max_x*d3_l->max_y;
194        carry=0;
195        if(!(*(d3_l->status+off)&AMORPH))
196        {
197         carry=(int)(my_info->dr_cc*(*(d3_l->extra+off)-*(d3_l->extra+offset))/2);
198         if(carry!=0)
199         {
200          *(d3_l->extra+offset)+=carry;
201          *(d3_l->extra+off)-=carry;
202         }
203        }
204       }
205      }
206      for(c=-1;c<=1;c++)
207      {
208       if(c!=0)
209       {
210        off=((i+c+d3_l->max_x)%d3_l->max_x)+j*d3_l->max_x+k*d3_l->max_x*d3_l->max_y;
211        carry=0;
212        if(!(*(d3_l->status+off)&AMORPH))
213        {
214         carry=(int)(my_info->dr_cc*(*(d3_l->extra+off)-*(d3_l->extra+offset))/2);
215         if(carry!=0)
216         {
217          *(d3_l->extra+offset)+=carry;
218          *(d3_l->extra+off)-=carry;
219         }
220        }
221       }
222      }
223     }
224    } /* for z */
225   } /* for y */
226  } /* for x */
227
228  } /* if step modulo diff_rate == 0 */
229
230  return 1;
231 }
232
233 int calc_pressure(d3_lattice *d3_l,int range)
234 {
235  int i,j,off;
236  double count,max=0;
237  int x,y,z;
238
239  for(x=0;x<d3_l->max_x;x++)
240  {
241   for(y=0;y<d3_l->max_y;y++)
242   {
243    for(z=0;z<d3_l->max_z;z++)
244    {
245     count=0;
246     for(i=-range;i<=range;i++)
247     {
248      for(j=-range;j<=range;j++)
249      {
250       if(i!=0 && j!=0)
251       {
252        off=((x+d3_l->max_x+i)%d3_l->max_x)+((y+d3_l->max_y+j)%d3_l->max_x)*d3_l->max_x+z*d3_l->max_x*d3_l->max_y;
253        if(*(d3_l->status+off)&AMORPH) count+=((double)*(d3_l->extra+off))/(i*i+j*j);
254       }
255      }
256     }
257     if(count>max) max=count;
258    }
259   }
260  }
261
262  for(x=0;x<d3_l->max_x;x++)
263  {
264   for(y=0;y<d3_l->max_y;y++)
265   {
266    for(z=0;z<d3_l->max_z;z++)
267    {
268     count=0;
269     for(i=-range;i<=range;i++)
270     {
271      for(j=-range;j<=range;j++)
272      {
273       if(i!=0 && j!=0)
274       {
275        off=((x+d3_l->max_x+i)%d3_l->max_x)+((y+d3_l->max_y+j)%d3_l->max_x)*d3_l->max_x+z*d3_l->max_x*d3_l->max_y;
276        if(*(d3_l->status+off)&AMORPH) count+=((double)*(d3_l->extra+off))/(i*i+j*j);
277       }
278      }
279     }
280     *(unsigned char *)(d3_l->v_ptr+x+y*d3_l->max_x+z*d3_l->max_x*d3_l->max_y)=(unsigned char)(count*255/max);
281    }
282   }
283  }
284
285  return 1;
286 }
287
288 int calc_max_extra(d3_lattice *d3_l)
289 {
290  int x,y,z;
291  int off,max=0;
292
293  for(x=0;x<d3_l->max_x;x++)
294  {
295   for(y=0;y<d3_l->max_y;y++)
296   {
297    for(z=0;z<d3_l->max_z;z++)
298    {
299     off=x+y*d3_l->max_x+z*d3_l->max_x*d3_l->max_y;
300     if(*(d3_l->extra+off)>max) max=*(d3_l->extra+off);
301    }
302   }
303  }
304
305  return max;
306 }
307
308 int write_bmp(d3_lattice *d3_l,int window,u32 x,u32 y,u32 z)
309 {
310  int fd,i,j,size=0;
311  int width=0,height=0;
312  char bmpfile[MAX_CHARS];
313  char buf[128];
314
315  if(window==1)
316  {
317   sprintf(bmpfile,"x-z_%d.bmp",y);
318   size=x*z;
319   width=x;
320   height=z;
321  }
322  if(window==2)
323  {
324   sprintf(bmpfile,"y-z_%d.bmp",x);
325   size=y*z;
326   width=y;
327   height=z;
328  }
329  if(window==3)
330  {
331   sprintf(bmpfile,"x-y_%d.bmp",z);
332   size=x*y;
333   width=x;
334   height=y;
335  }
336
337  if((fd=open(bmpfile,O_WRONLY|O_CREAT))<0)
338  {
339   puts("cannot open bmp file");
340   return -1;
341  }
342
343  /* bmpheader */
344  buf[0]='B'; /* std header start */
345  buf[1]='M';
346  buf[2]=(size+0x36)&0xff; /* file size */
347  buf[3]=(size+0x36)>>8;
348  memset(buf+4,0,6);
349  buf[10]=0x36; /* offset to data */
350  memset(buf+11,0,3);
351  buf[14]=0x28; /* length of bmp info header */
352  memset(buf+15,0,3);
353  buf[18]=width&0xff; /* width and height */
354  buf[19]=width>>8;
355  memset(buf+20,0,2);
356  buf[22]=height&0xff;
357  buf[23]=height>>8;
358  memset(buf+24,0,2);
359  buf[26]=1; /* # planes -> 1 */
360  buf[27]=0;
361  buf[28]=24; /* bits per pixel -> 2^24 (true color) */
362  buf[29]=0;
363  memset(buf+30,0,4); /* compression -> none */
364  buf[34]=size&0xff; /* data size */
365  buf[35]=size>>8;
366  memset(buf+36,0,2);
367  buf[38]=0x12; /* res: pixel/meter */
368  buf[39]=0x0b;
369  memset(buf+40,0,2);
370  buf[42]=0x12;
371  buf[43]=0x0b;
372  memset(buf+44,0,2); 
373  memset(buf+46,0,8); /* no colors, no important colors */
374
375  if(write(fd,buf,54)<54)
376  {
377   puts("failed writing bmp header");
378   return -1;
379  }
380  if(window==1)
381  {
382   for(i=0;i<d3_l->max_x;i++)
383   {
384    for(j=0;j<d3_l->max_z;j++)
385    {
386     if(*(d3_l->status+i+y*d3_l->max_x+j*d3_l->max_x*d3_l->max_y)&RED) memset(buf,0xff,3);
387     else memset(buf,0,3);
388     if(write(fd,buf,3)<3)
389     {
390      puts("failed writing rgb values to bmp file");
391      return-1;
392     }
393    }
394   }
395  }
396  if(window==2)
397  {
398   for(i=0;i<d3_l->max_y;i++)
399   {
400    for(j=0;j<d3_l->max_z;j++)
401    {
402     if(*(d3_l->status+x+i*d3_l->max_x+j*d3_l->max_x*d3_l->max_y)&RED) memset(buf,0xff,3);
403     else memset(buf,0,3);
404     if(write(fd,buf,3)<3)
405     {
406      puts("failed writing rgb values to bmp file");
407      return-1;
408     }
409    }
410   }
411  }
412  if(window==3)
413  {
414   for(i=0;i<d3_l->max_x;i++)
415   {
416    for(j=0;j<d3_l->max_y;j++)
417    {
418     if(*(d3_l->status+i+j*d3_l->max_x+z*d3_l->max_x*d3_l->max_y)&RED) memset(buf,0xff,3);
419     else memset(buf,0,3);
420     if(write(fd,buf,3)<3)
421     {
422      puts("failed writing rgb values to bmp file");
423      return -1;
424     }
425    }
426   }
427  }
428  close(fd);
429
430  return 1;
431 }
432
433 int save_to_file(char *sf,d3_lattice *d3_l,info *my_inf)
434 {
435  int sf_fd,c;
436
437  if((sf_fd=open(sf,O_WRONLY|O_CREAT))<0)
438  {
439   puts("cannot open save file");
440   return -1;
441  }
442  if(write(sf_fd,d3_l,sizeof(d3_lattice))<sizeof(d3_lattice))
443  {
444   puts("failed saving d3 lattice struct");
445   return -1;
446  }
447  if(write(sf_fd,my_inf,sizeof(info))<sizeof(info))
448  {
449   puts("failed saving info struct");
450   return-1;
451  }
452  c=d3_l->max_x*d3_l->max_y*d3_l->max_z;
453  if(write(sf_fd,d3_l->status,c*sizeof(unsigned char))<c*sizeof(unsigned char))
454  {
455   puts("failed saving status of d3 lattice sites");
456   return -1;
457  }
458  if(write(sf_fd,d3_l->extra,c*sizeof(int))<c*sizeof(int))
459  {
460   puts("failed saving sites concentration");
461   return -1;
462  }
463  close(sf_fd);
464
465  return 1;
466 }
467
468 int load_from_file(char *lf,d3_lattice *d3_l,info *my_inf)
469 {
470  int lf_fd,c;
471
472  if((lf_fd=open(lf,O_RDONLY))<0)
473  {
474   puts("cannot open load file");
475   return -1;
476  }
477  if(read(lf_fd,d3_l,sizeof(d3_lattice))<sizeof(d3_lattice))
478  {
479   puts("failed reading d3 lattice struct");
480   return -1;
481  }
482  if(read(lf_fd,my_inf,sizeof(info))<sizeof(info))
483  {
484   puts("failed reading info struct");
485   return-1;
486  }
487  c=d3_l->max_x*d3_l->max_y*d3_l->max_z;
488  if((d3_l->status=(unsigned char*)malloc(c*sizeof(unsigned char)))==NULL)
489  {
490   puts("cannot allocate status buffer");
491   return -1;
492  }
493  if((d3_l->extra=(int *)malloc(c*sizeof(int)))==NULL)
494  {
495   puts("cannot allocate concentration buffer");
496   return -1;
497  }
498  if(read(lf_fd,d3_l->status,c*sizeof(unsigned char))<c*sizeof(unsigned char))
499  {
500   puts("failed reading status of d3 lattice sites");
501   return -1;
502  }
503  if(read(lf_fd,d3_l->extra,c*sizeof(int))<c*sizeof(int))
504  {
505   puts("failed reading sites concentration");
506   return -1;
507  }
508  close(lf_fd);
509
510  return 1;
511 }
512
513 int convert_file(char *cf,d3_lattice *d3_l)
514 {
515  int x,y,z;
516  int c_fd;
517
518  if((c_fd=open(cf,O_WRONLY|O_CREAT))<0)
519  {
520   puts("cannot open convert file");
521   return -1;
522  }
523  dprintf(c_fd,"# created by nlsop (gnuplot format)\n");
524  for(x=0;x<d3_l->max_x;x++)
525  {
526   for(y=0;y<d3_l->max_y;y++)
527   {
528    for(z=0;z<d3_l->max_z;z++)
529    {
530     if(*(d3_l->status+x+y*d3_l->max_x+z*d3_l->max_x*d3_l->max_y)&AMORPH) dprintf(c_fd,"%d %d %d\n",x,y,z);
531    }
532   }
533  }
534  close(c_fd);
535
536  return 1;
537 }
538
539 int get_c_ratio(double *c_ratio,char *pfile,info *my_info,d3_lattice *d3_l)
540 {
541  double all,a,b,d;
542  int i,k;
543  int p_fd;
544  unsigned char buf[32];
545  char *p;
546  unsigned char c;
547
548  if((p_fd=open(pfile,O_RDONLY))<0)
549  {
550   puts("cannot open profile file");
551   return -1;
552  }
553  k=1;
554  d=0;
555  all=0;
556  while(k)
557  {
558   for(i=0;i<32;i++)
559   {
560    k=read(p_fd,&c,1);
561    buf[i]=c;
562    if(c=='\n') break;
563   }
564   if(k)
565   {
566    p=strtok(buf," ");
567    a=atof(p)/10; /* nm */
568    p=strtok(NULL," ");
569    b=atof(p);
570    if(a>my_info->b_cd*CELL_LENGTH && a<(my_info->b_cd+d3_l->max_z)*CELL_LENGTH) d+=b;
571    all+=b;
572   }
573  }
574  *c_ratio=d/all;
575  close(p_fd);
576
577  return 1;
578 }
579
580 int main(int argc,char **argv)
581 {
582  u32 x,y,z,x_c,y_c,z_c;
583  int i,quit,escape,switchmode,nowait,bmp;
584  int refresh,resave;
585  char s_file[MAX_CHARS];
586  char s_file_tmp[MAX_CHARS];
587  char l_file[MAX_CHARS];
588  char c_file[MAX_CHARS];
589  char p_file[MAX_CHARS];
590  char convert;
591  char r_file[MAX_CHARS];
592 #ifdef USE_DFB_API
593  char xyz_txt[MAX_TXT];
594  char status_txt[MAX_TXT];
595  char conc_txt[MAX_TXT];
596  char steps_txt[MAX_TXT];
597  char cc_txt[MAX_TXT];
598  char a_txt[MAX_TXT];
599  char s_txt[MAX_TXT];
600  char ap_txt[MAX_TXT];
601  char el_txt[MAX_TXT];
602  char cd_txt[MAX_TXT];
603  char r_txt[MAX_TXT];
604  char cp_txt[MAX_TXT];
605  char zdiff_txt[MAX_TXT];
606  char diff_txt[MAX_TXT];
607  char dr_ac_txt[MAX_TXT];
608  char dr_cc_txt[MAX_TXT];
609  char dose_txt[MAX_TXT];
610  char mode_txt[MAX_TXT];
611  char *arg_v[MAX_ARGV];
612 #endif
613  d3_lattice d3_l;
614  info my_info;
615  unsigned char mode;
616  double c_ratio;
617 #ifdef USE_DFB_API
618  int max_extra;
619 #endif
620
621  d3_l.max_x=X;
622  d3_l.max_y=Y;
623  d3_l.max_z=Z;
624  my_info.steps=STEPS;
625  my_info.range=RANGE;
626  refresh=REFRESH;
627  resave=RESAVE;
628  my_info.z_diff=0;
629  my_info.a_el=A_EL;
630  my_info.b_el=B_EL;
631  my_info.a_cd=A_CD;
632  my_info.b_cd=B_CD;
633  my_info.a_ap=A_AP;
634  my_info.b_ap=B_AP;
635  my_info.a_cp=A_CP;
636  my_info.cc=CC;
637  my_info.dr_ac=DR_AC;
638  my_info.dr_cc=DR_CC;
639  my_info.diff_rate=DIFF_RATE;
640  nowait=0;
641  quit=0;
642  escape=0;
643  switchmode=0;
644  strcpy(s_file,"");
645  strcpy(l_file,"");
646  strcpy(c_file,"");
647  strcpy(p_file,IMP_PROFILE);
648  convert=0;
649  strcpy(r_file,"");
650  mode=0;
651
652  amorph_count=0;
653  cryst_count=0;
654
655  for(i=1;i<argc;i++)
656  {
657   if(argv[i][0]=='-')
658   {
659    switch(argv[i][1])
660    {
661     case 'h':
662      usage();
663      return -1;
664     case 'n':
665      nowait=1;
666      break;
667     case 'a':
668      my_info.a_el=atof(argv[++i]);
669      break;
670     case 'b':
671      my_info.b_el=atof(argv[++i]);
672      break;
673     case 'x':
674      d3_l.max_x=atoi(argv[++i]);
675      break;
676     case 'y':
677      d3_l.max_y=atoi(argv[++i]);
678      break;
679     case 'z':
680      d3_l.max_z=atoi(argv[++i]);
681      break;
682     case 'Z':
683      my_info.z_diff=1;
684      break;
685     case 's':
686      my_info.steps=atoi(argv[++i]);
687      break;
688     case 'd':
689      refresh=atoi(argv[++i]);
690      break;
691     case 'r':
692      my_info.range=atoi(argv[++i]);
693      break;
694     case 'f':
695      my_info.a_ap=atof(argv[++i]);
696      break;
697     case 'p':
698      my_info.b_ap=atof(argv[++i]);
699      break;
700     case 'F':
701      my_info.a_cp=atof(argv[++i]);
702      break;
703     case 'A':
704      my_info.a_cd=atof(argv[++i]);
705      break;
706     case 'B':
707      my_info.b_cd=atof(argv[++i]);
708      break;
709     case 'W':
710      resave=atoi(argv[++i]);
711      break;
712     case 'C':
713      strcpy(l_file,argv[++i]);
714      if(i<argc-1) if(argv[i+1][0]!='-') strcpy(c_file,argv[++i]);
715      convert=1;
716      break;
717     case 'D':
718      my_info.dr_ac=atof(argv[++i]);
719      break;
720     case 'c':
721      my_info.dr_cc=atof(argv[++i]);
722      break;
723     case 'e':
724      my_info.diff_rate=atoi(argv[++i]);
725      break;
726     case 'L':
727      strcpy(l_file,argv[++i]);
728      break;
729     case 'S':
730      strcpy(s_file,argv[++i]);
731      break;
732     case 'R':
733      strcpy(r_file,argv[++i]);
734      break;
735     case 'P':
736      strcpy(p_file,argv[++i]);
737      break;
738     default:
739      usage();
740      return -1;
741    }
742   } else usage();
743  }
744
745  x=d3_l.max_x/2-1;
746  y=d3_l.max_y/2-1;
747  z=d3_l.max_z/2-1;
748
749 #ifdef NODFB
750  if(!strcmp(s_file,""))
751  {
752   puts("NODFB defined, run with -S option");
753   return -1;
754  }
755 #endif
756
757  if(!strcmp(r_file,"")) rand_init(NULL);
758  else rand_init(r_file);
759
760  if(!strcmp(l_file,""))
761  {
762   i=d3_l.max_x*d3_l.max_y*d3_l.max_z;
763 #ifdef USE_DFB_API
764   d3_lattice_init(&argc,argv,&d3_l);
765 #endif
766   if((d3_l.status=(unsigned char *)malloc(i*sizeof(unsigned char)))==NULL)
767   {
768    puts("failed allocating status buffer");
769    return -1;
770   }
771   memset(d3_l.status,0,i*sizeof(unsigned char));
772   if((d3_l.extra=(int *)malloc(i*sizeof(int)))==NULL)
773   {
774    puts("failed allocating concentration buffer");
775    return -1;
776   }
777   memset(d3_l.extra,0,i*sizeof(int));
778  } else
779  {
780   load_from_file(l_file,&d3_l,&my_info);
781   if(convert) 
782   {   
783    if(!strcmp(c_file,"")) sprintf(c_file,"%s_gnuplot",l_file);
784    printf("converting file %s to %s\n",l_file,c_file);
785    convert_file(c_file,&d3_l);
786    puts("done");
787    return 1;
788   } 
789 #ifdef USE_DFB_API
790     else d3_lattice_init(&argc,argv,&d3_l);
791 #endif
792  }
793
794 #ifdef USE_DFB_API
795  d3_event_init(&d3_l);
796 #endif
797
798 #ifdef USE_DFB_API
799  strcpy(a_txt,"args:");
800  sprintf(s_txt,"steps: %d",my_info.steps);
801  sprintf(dose_txt,"dose: %.2fe+17 C/cm²",my_info.steps*1.0/(d3_l.max_x*d3_l.max_y*CELL_LENGTH*CELL_LENGTH*1000));
802  sprintf(r_txt,"pressure range: %d",my_info.range);
803  sprintf(ap_txt,"a_ap: %.3f  b_ap: %.3f",my_info.a_ap,my_info.b_ap);
804  sprintf(el_txt,"a_el: %.3f  b_el: %.3f",my_info.a_el,my_info.b_el);
805  sprintf(cd_txt,"a_cd: %.3f  b_cd: %.3f",my_info.a_cd,my_info.b_cd);
806  sprintf(cp_txt,"a_cp: %.4f",my_info.a_cp);
807  sprintf(dr_ac_txt,"a/c diffusion rate: %.4f",my_info.dr_ac);
808  sprintf(dr_cc_txt,"c/c diffusion rate: %.4f",my_info.dr_cc);
809  sprintf(zdiff_txt,"diffusion in z direction: %c",my_info.z_diff?'y':'n');
810  sprintf(diff_txt,"diffusion every %d steps",my_info.diff_rate);
811  strcpy(mode_txt,"view: a/c mode");
812  arg_v[1]=xyz_txt;
813  arg_v[2]=NULL;
814  arg_v[3]=status_txt;
815  arg_v[4]=conc_txt;
816  arg_v[5]=NULL;
817  arg_v[6]=mode_txt;
818  arg_v[7]=NULL;
819  arg_v[8]=steps_txt;
820  arg_v[9]=cc_txt;
821  arg_v[10]=NULL;
822  arg_v[11]=diff_txt;
823  arg_v[12]=zdiff_txt;
824  arg_v[13]=NULL;
825  arg_v[14]=a_txt;
826  arg_v[15]=s_txt;
827  arg_v[16]=dose_txt;
828  arg_v[17]=NULL;
829  arg_v[18]=r_txt;
830  arg_v[19]=ap_txt;
831  arg_v[20]=el_txt;
832  arg_v[21]=cd_txt;
833  arg_v[22]=cp_txt;
834  arg_v[23]=NULL;
835  arg_v[24]=dr_ac_txt;
836  arg_v[25]=dr_cc_txt;
837 #endif
838
839  if(!strcmp(l_file,""))
840  {
841   if(get_c_ratio(&c_ratio,p_file,&my_info,&d3_l)!=1)
842   {
843    puts("failed calculating ratio");
844    return -1;
845   }
846   i=0;
847   while((i<my_info.steps) && (quit==0) && (escape==0))
848   {
849    x_c=get_rand(d3_l.max_x);
850    y_c=get_rand(d3_l.max_y);
851    z_c=get_rand_lgp(d3_l.max_z,my_info.a_el,my_info.b_el);
852    distrib_c(&d3_l,&my_info,i,c_ratio);
853    process_cell(&d3_l,x_c,y_c,z_c,&my_info);
854 #ifdef USE_DFB_API
855    if(i%refresh==0)
856    {
857     sprintf(xyz_txt,"x: %d  y: %d  z: %d",x+1,y+1,z+1);
858     sprintf(status_txt,"status: %c",(*(d3_l.status+x+y*d3_l.max_x+z*d3_l.max_x*d3_l.max_y)&AMORPH)?'a':'c');
859     sprintf(conc_txt,"conc: %d",*(d3_l.extra+x+y*d3_l.max_x+z*d3_l.max_x*d3_l.max_y));
860     sprintf(steps_txt,"step: %d  a_count: %d",i,amorph_count);
861     sprintf(cc_txt,"total c: %d  c_count: %d",my_info.cc,cryst_count);
862     d3_lattice_draw(&d3_l,x,y,z,25,arg_v,mode,0);
863    }
864 #endif
865    if(i%resave==0 && strcmp(s_file,"") && resave!=0 && i!=0)
866    {
867     sprintf(s_file_tmp,"%s_%d_of_%d",s_file,i,my_info.steps);
868     save_to_file(s_file_tmp,&d3_l,&my_info);
869 #ifdef NODFB
870     printf("saved %s\n",s_file_tmp);
871 #endif
872    }
873    i++;
874   }
875  }
876
877  if(strcmp(s_file,""))
878  {
879    printf("saved %s\n",s_file);
880    save_to_file(s_file,&d3_l,&my_info);
881  }
882
883 #ifdef USE_DFB_API
884  /* allocating buffer for pressure values */
885  if((d3_l.v_ptr=malloc(d3_l.max_x*d3_l.max_y*d3_l.max_z*sizeof(unsigned char)))==NULL)
886  {
887   puts("cannot allocate buffer for pressure values");
888   return -1;
889  }
890  /* calc values */
891  calc_pressure(&d3_l,my_info.range);
892  max_extra=calc_max_extra(&d3_l);
893
894  while((quit==0) && (escape==0) && (nowait==0))
895  {
896   /* bahh! */
897   if(switchmode==0) mode=0;
898   if(switchmode==1) mode=1;
899   if(switchmode==2) mode=2;
900   /* end of bahh! */
901   sprintf(xyz_txt,"x: %d  y: %d  z: %d",x+1,y+1,z+1);
902   sprintf(status_txt,"status: %c",(*(d3_l.status+x+y*d3_l.max_x+z*d3_l.max_x*d3_l.max_y)&AMORPH)?'a':'c');
903   sprintf(conc_txt,"conc: %d",*(d3_l.extra+x+y*d3_l.max_x+z*d3_l.max_x*d3_l.max_y));
904   sprintf(steps_txt,"step: end");
905   sprintf(cc_txt,"total c: %d",my_info.cc);
906   if(switchmode==0) strcpy(mode_txt,"view: a/c mode");
907   if(switchmode==1) strcpy(mode_txt,"view: c conc mode");
908   if(switchmode==2) strcpy(mode_txt,"view: a pressure mode");
909   d3_lattice_draw(&d3_l,x,y,z,25,arg_v,mode,max_extra);
910   bmp=0;
911   scan_event(&d3_l,&x,&y,&z,&quit,&escape,&switchmode,&bmp);
912   if(bmp) write_bmp(&d3_l,bmp,x,y,z);
913
914  }
915
916  d3_lattice_release(&d3_l);
917 #endif
918
919  return 1;
920 }