ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/group/trunk/xyz2pov/src/xyz2pov.c
Revision: 868
Committed: Wed Nov 19 16:26:33 2003 UTC (20 years, 7 months ago) by gezelter
Content type: text/plain
File size: 34116 byte(s)
Log Message:
Fixes for rotation

File Contents

# User Rev Content
1 mmeineke 60 #include <stdio.h>
2     #include <stdlib.h>
3     #include <string.h>
4     #include <math.h>
5     #include <unistd.h>
6     #include <sys/types.h>
7     #include <sys/stat.h>
8    
9     #include "frameCount.h"
10     #include "atom_parser.h"
11     #include "pov_writer.h"
12    
13    
14     #define POV_DIR "./pov"
15    
16 mmeineke 508
17 mmeineke 60 struct linked_xyz{
18     struct coords *r;
19 gezelter 625 double Hmat[3][3];
20 mmeineke 60 struct linked_xyz *next;
21     };
22    
23     char *program_name; /*the name of the program */
24     int draw_bonds = 0; /* boolean to draw bonds or not */
25     int draw_hydrogens = 0; /*boolean to draw hydrogens */
26     int draw_atoms = 0; /*boolean to draw atoms */
27 gezelter 864 int draw_vectors = 0; /*boolean to draw vectors */
28 mmeineke 515 int draw_box = 0; // boolean to draw the periodic Box
29 mmeineke 508 int regenerateBonds = 0; // boolean to regenearate bonds each frame
30 mmeineke 60
31     void usage(void);
32 gezelter 864 int count_tokens(char *line, char *delimiters);
33 mmeineke 60
34     int main(argc, argv)
35     int argc;
36     char *argv[];
37     {
38    
39    
40     struct coords *out_coords;
41    
42 gezelter 625 int i,j,k; /* loop counters */
43 mmeineke 60 mode_t dir_mode = S_IRWXU;
44    
45     int generate_header = 0; /* boolean for generating the pov ray header */
46     double big_x = 0;
47     double big_y = 0; /* lets me know the biggest x y and z */
48     double big_z = 0;
49     double small_x = 0;
50     double small_y = 0; /* lets me know the smallest x, y, z */
51     double small_z = 0;
52 gezelter 864 int extremaSet = 0;
53 mmeineke 60 double rsqr; /* the square of the diagonal */
54     double diagonal; /* the diagonal length of the sim box */
55    
56     unsigned int n_atoms; /*the number of atoms in each time step */
57 mmeineke 649 char read_buffer[2000]; /*the line buffer for reading */
58 mmeineke 60 char *eof_test; /*ptr to see when we reach the end of the file */
59     char *foo; /*the pointer to the current string token */
60     FILE *in_file; /* the input file */
61     FILE *out_file; /*the output file */
62     char *out_prefix = NULL; /*the prefix of the output file */
63     int have_prefix = 0;
64     char out_name[500]; /*the output name */
65     char out_format[1000];
66     char *in_name = NULL; /*the name of the input file */
67     unsigned int n_out = 0; /*keeps track of which output file is being written*/
68     int done;
69     char current_flag;
70     int nFrames;
71     int nZeroes;
72 gezelter 864 int nTokens;
73 mmeineke 60 double count;
74    
75 mmeineke 515 int startFrame = 1;
76     int endFrame;
77     int span = 0;
78     int currentCount = 0;
79     int haveEnd = 0;
80    
81 mmeineke 60 struct linked_xyz *current_frame;
82     struct linked_xyz *temp_frame;
83    
84     unsigned int n_interpolate = 0; /* number of frames to interpolate */
85     double dx, dy, dz; /* temp variables for interpolating distances */
86 gezelter 625 double dm[3][3];
87 mmeineke 515
88 mmeineke 60
89     char pov_dir[500]; /* the pov_dir */
90    
91     program_name = argv[0]; /*save the program name in case we need it*/
92    
93 mmeineke 508
94 mmeineke 60 for( i = 1; i < argc; i++){
95    
96     if(argv[i][0] =='-'){
97    
98     // parse the option
99    
100     if(argv[i][1] == '-' ){
101    
102     // parse long word options
103    
104     fprintf( stderr,
105     "Invalid option \"%s\"\n", argv[i] );
106     usage();
107    
108     }
109    
110     else{
111    
112     // parse single character options
113    
114     done = 0;
115     j = 1;
116     current_flag = argv[i][j];
117     while( (current_flag != '\0') && (!done) ){
118    
119     switch(current_flag){
120    
121     case 'o':
122     // -o <prefix> => the output prefix.
123    
124     i++;
125     out_prefix = argv[i];
126     have_prefix = 1;
127     done = 1;
128     break;
129    
130     case 'i':
131     // -i <#> => the number to interpolate
132    
133     i++;
134     n_interpolate = atoi( argv[i] );
135     done = 1;
136     break;
137    
138 mmeineke 515 case 'f':
139     // -f <#> => frame to render
140    
141     i++;
142     startFrame = atoi( argv[i] );
143     endFrame = startFrame;
144     haveEnd = 1;
145     done = 1;
146     break;
147    
148     case 's':
149     // -s <#> => frame to start
150    
151     i++;
152     startFrame = atoi( argv[i] );
153     done = 1;
154     break;
155    
156     case 'e':
157     // -e <#> => frame to end
158    
159     i++;
160     endFrame = atoi( argv[i] );
161     haveEnd = 1;
162     done = 1;
163     break;
164    
165 mmeineke 60 case 'H':
166     // -h => generate a pov-ray Header
167    
168     generate_header = 1;
169     break;
170    
171     case 'h':
172     // -h => draw Hydrogens
173    
174     draw_hydrogens = 1;
175     break;
176    
177     case 'b':
178     // -b => draw bonds
179    
180     draw_bonds = 1;
181     break;
182    
183     case 'a':
184     // -a => draw the atoms
185    
186     draw_atoms = 1;
187     break;
188    
189 gezelter 864 case 'v':
190     // -v => draw the vectors
191    
192     draw_vectors = 1;
193     break;
194    
195 mmeineke 508 case 'r':
196     // -r => regenerate bonds
197    
198     regenerateBonds = 1;
199     break;
200    
201 mmeineke 515 case 'p':
202     // -r => draw periodic box
203    
204     draw_box = 1;
205     break;
206    
207 mmeineke 60 default:
208    
209     (void)fprintf(stderr, "Bad option \"-%c\"\n", current_flag);
210     usage();
211     }
212     j++;
213     current_flag = argv[i][j];
214     }
215     }
216     }
217    
218     else{
219    
220     if( in_name != NULL ){
221     fprintf( stderr,
222     "Error at \"%s\", program does not currently support\n"
223     "more than one input file.\n"
224     "\n",
225     argv[i]);
226     usage();
227     }
228    
229     in_name = argv[i];
230     }
231     }
232    
233     if(in_name == NULL){
234     usage();
235     }
236    
237    
238    
239     in_file = fopen(in_name, "r");
240     if(in_file == NULL){
241     printf("Cannot open file: %s\n", in_name);
242     exit(8);
243     }
244    
245    
246    
247     if(access(POV_DIR, F_OK)){
248     /*create the pov directory*/
249     mkdir(POV_DIR, dir_mode);
250     }
251     strcpy(pov_dir, POV_DIR); strcat(pov_dir, "/");
252    
253    
254     // initialize atom type parser
255    
256     initializeParser();
257 mmeineke 508 initBondList();
258 mmeineke 60
259     // count the number of frames
260    
261     printf( "Counting the number of frames..." );
262     fflush(stdout);
263    
264     nFrames = frameCount( in_name );
265    
266     printf( "done.\n"
267     "%d frames found\n",
268     nFrames);
269     fflush(stdout);
270    
271     // create the output string
272    
273     nZeroes = 1;
274     count = (double)( nFrames * (n_interpolate+1) );
275     while( count >= 10.0 ){
276     count /= 10.0;
277     nZeroes++;
278     }
279    
280     if(!have_prefix){
281     out_prefix = strtok(in_name, ".");
282     }
283    
284 mmeineke 649 sprintf( out_format, "%s%s-%%0%dd.pov", pov_dir, out_prefix, nZeroes );
285 mmeineke 60
286     // start reading the first frame
287    
288     eof_test = fgets(read_buffer, sizeof(read_buffer), in_file);
289    
290     current_frame = (struct linked_xyz *)malloc(sizeof(struct linked_xyz));
291     current_frame->next = NULL;
292    
293    
294 mmeineke 515 if( haveEnd ) span = endFrame - startFrame;
295     done = 0;
296     if( span < 0 ) done = 1;
297     while( (eof_test != NULL) && !done ){
298 mmeineke 60
299     (void)sscanf(read_buffer, "%d", &n_atoms);
300     current_frame->r =
301     (struct coords *)calloc(n_atoms, sizeof(struct coords));
302    
303     /*read and toss the comment line */
304    
305     eof_test = fgets(read_buffer, sizeof(read_buffer), in_file);
306     if(eof_test == NULL){
307     printf("error in reading file\n");
308     exit(8);
309     }
310 mmeineke 515
311     // unless we need to get the box size
312    
313     if( draw_box ){
314     foo = strtok(read_buffer, " ,;\t");
315     if(foo == NULL){
316 mmeineke 649 printf("error in reading file time\n");
317 mmeineke 515 exit(8);
318     }
319    
320     foo = strtok(NULL, " ,;\t");
321     if(foo == NULL){
322 mmeineke 649 printf("error in reading file h00\n");
323 mmeineke 515 exit(8);
324     }
325 gezelter 625 current_frame->Hmat[0][0] = atof( foo );
326 mmeineke 60
327 mmeineke 515 foo = strtok(NULL, " ,;\t");
328     if(foo == NULL){
329 mmeineke 649 printf("error in reading file h10\n");
330 mmeineke 515 exit(8);
331     }
332 gezelter 625 current_frame->Hmat[1][0] = atof( foo );
333 mmeineke 515
334     foo = strtok(NULL, " ,;\t");
335     if(foo == NULL){
336 mmeineke 649 printf("error in reading file h20\n");
337 mmeineke 515 exit(8);
338     }
339 gezelter 625 current_frame->Hmat[2][0] = atof( foo );
340    
341     foo = strtok(NULL, " ,;\t");
342     if(foo == NULL){
343 mmeineke 649 printf("error in reading file h01\n");
344 gezelter 625 exit(8);
345     }
346     current_frame->Hmat[0][1] = atof( foo );
347    
348     foo = strtok(NULL, " ,;\t");
349     if(foo == NULL){
350 mmeineke 649 printf("error in reading file h11\n");
351 gezelter 625 exit(8);
352     }
353     current_frame->Hmat[1][1] = atof( foo );
354    
355     foo = strtok(NULL, " ,;\t");
356     if(foo == NULL){
357 mmeineke 649 printf("error in reading file h21\n");
358 gezelter 625 exit(8);
359     }
360     current_frame->Hmat[2][1] = atof( foo );
361    
362     foo = strtok(NULL, " ,;\t");
363     if(foo == NULL){
364 mmeineke 649 printf("error in reading file h02\n");
365 gezelter 625 exit(8);
366     }
367     current_frame->Hmat[0][2] = atof( foo );
368    
369     foo = strtok(NULL, " ,;\t");
370     if(foo == NULL){
371 mmeineke 649 printf("error in reading file h12\n");
372 gezelter 625 exit(8);
373     }
374     current_frame->Hmat[1][2] = atof( foo );
375    
376     foo = strtok(NULL, " ,;\t");
377     if(foo == NULL){
378 mmeineke 649 printf("error in reading file h22\n");
379 gezelter 625 exit(8);
380     }
381     current_frame->Hmat[2][2] = atof( foo );
382    
383 mmeineke 515 }
384    
385 mmeineke 60 for( i=0; i < n_atoms; i++){
386    
387     eof_test = fgets(read_buffer, sizeof(read_buffer), in_file);
388     if(eof_test == NULL){
389 mmeineke 649 printf("error in reading file line at atom %d\n", i);
390 mmeineke 60 exit(8);
391     }
392    
393 gezelter 864 nTokens = count_tokens(read_buffer, " ,;\t");
394    
395     if (nTokens < 4) {
396     printf("Not enough tokens while parsing file at atom %d\n", i);
397     exit(8);
398     }
399    
400 mmeineke 60 foo = strtok(read_buffer, " ,;\t");
401     (void)strcpy(current_frame->r[i].name, foo); /*copy the atom name */
402    
403     foo = strtok(NULL, " ,;\t");
404     (void)sscanf(foo, "%lf",&current_frame->r[i].x);
405     foo = strtok(NULL, " ,;\t");
406     (void)sscanf(foo, "%lf", &current_frame->r[i].y);
407 gezelter 864 foo = strtok(NULL, " ,;\t");
408     (void)sscanf(foo, "%lf", &current_frame->r[i].z);
409 mmeineke 60
410 gezelter 864 if (extremaSet) {
411     if(current_frame->r[i].x > big_x) big_x = current_frame->r[i].x;
412     if(current_frame->r[i].x < small_x) small_x = current_frame->r[i].x;
413    
414     if(current_frame->r[i].y > big_y) big_y = current_frame->r[i].y;
415     if(current_frame->r[i].y < small_y) small_y = current_frame->r[i].y;
416    
417     if(current_frame->r[i].z > big_z) big_z = current_frame->r[i].z;
418     if(current_frame->r[i].z < small_z) small_z = current_frame->r[i].z;
419     } else {
420     big_x = current_frame->r[i].x;
421     small_x = current_frame->r[i].x;
422    
423     big_y = current_frame->r[i].y;
424     small_y = current_frame->r[i].y;
425    
426     big_z = current_frame->r[i].z;
427     small_z = current_frame->r[i].z;
428 mmeineke 60
429 gezelter 864 extremaSet = 1;
430    
431 mmeineke 60 }
432    
433 gezelter 864 if (nTokens == 5 || nTokens > 7) {
434     foo = strtok(NULL, " ,;\t");
435     (void)sscanf(foo, "%lf", &current_frame->r[i].charge);
436     current_frame->r[i].hasCharge = 1;
437     } else {
438     current_frame->r[i].hasCharge = 0;
439     }
440    
441    
442     if (nTokens >= 7) {
443     foo = strtok(NULL, " ,;\t");
444     (void)sscanf(foo, "%lf", &current_frame->r[i].ux);
445     foo = strtok(NULL, " ,;\t");
446     (void)sscanf(foo, "%lf", &current_frame->r[i].uy);
447     foo = strtok(NULL, " ,;\t");
448     (void)sscanf(foo, "%lf", &current_frame->r[i].uz);
449     current_frame->r[i].hasVector = 1;
450     } else {
451     current_frame->r[i].hasVector = 0;
452     }
453 mmeineke 60 }
454 mmeineke 515 currentCount++;
455 mmeineke 60
456 mmeineke 515
457     if( currentCount >= startFrame ){
458     if(n_interpolate && current_frame->next != NULL){
459 mmeineke 60
460 mmeineke 515 temp_frame = current_frame->next;
461 mmeineke 60
462 mmeineke 515 for(i = 0; i < n_interpolate; i++){
463    
464     /* open the new output file */
465    
466     sprintf(out_name, out_format, currentCount );
467     out_file = fopen(out_name, "w");
468     currentCount++;
469     if(out_file == NULL){
470     printf("error opening output file: %s\n", out_name);
471     exit(8);
472     }
473     (void)fprintf(out_file,
474     "// The following script was automatically generated by:\n"
475     "// xyz2pov Copyright 2001 by MATTHEW A. MEINEKE\n"
476     "\n"
477     "#include \"pov_header.pov\"\n"
478     "\n");
479     if( draw_box ){
480 gezelter 625
481     for (j = 0; j < 3; j++) {
482     for (k = 0; k < 3; k++) {
483     dm[j][k] = current_frame->Hmat[j][k] - temp_frame->Hmat[j][k];
484     dm[j][k] /= (double)(n_interpolate + 1);
485     }
486     }
487    
488 mmeineke 515 fprintf( out_file,
489 gezelter 625 "makePeriodicBox( %lf, %lf, %lf, %lf, %lf, %lf, %lf, %lf, %lf)\n"
490 mmeineke 515 "\n",
491 gezelter 625 temp_frame->Hmat[0][0] + dm[0][0] * (i+1),
492 mmeineke 629 temp_frame->Hmat[2][0] + dm[2][0] * (i+1),
493 gezelter 625 temp_frame->Hmat[1][0] + dm[1][0] * (i+1),
494     temp_frame->Hmat[0][1] + dm[0][1] * (i+1),
495 mmeineke 629 temp_frame->Hmat[2][1] + dm[2][1] * (i+1),
496 gezelter 625 temp_frame->Hmat[1][1] + dm[1][1] * (i+1),
497     temp_frame->Hmat[0][2] + dm[0][2] * (i+1),
498 mmeineke 629 temp_frame->Hmat[2][2] + dm[2][2] * (i+1),
499     temp_frame->Hmat[1][2] + dm[1][2] * (i+1) );
500 mmeineke 515 }
501    
502    
503     out_coords =
504     (struct coords *)calloc(n_atoms, sizeof(struct coords));
505    
506     for(j=0; j < n_atoms; j++){
507     dx = current_frame->r[j].x - temp_frame->r[j].x;
508     dy = current_frame->r[j].y - temp_frame->r[j].y;
509     dz = current_frame->r[j].z - temp_frame->r[j].z;
510    
511     dx /= (double)(n_interpolate + 1);
512     dy /= (double)(n_interpolate + 1);
513     dz /= (double)(n_interpolate + 1);
514    
515     strcpy(out_coords[j].name, temp_frame->r[j].name);
516     out_coords[j].x = temp_frame->r[j].x + dx * (i+1);
517     out_coords[j].y = temp_frame->r[j].y + dy * (i+1);
518     out_coords[j].z = temp_frame->r[j].z + dz * (i+1);
519 gezelter 864
520     if (current_frame->r[j].hasVector) {
521     dx = current_frame->r[j].ux - temp_frame->r[j].ux;
522     dy = current_frame->r[j].uy - temp_frame->r[j].uy;
523     dz = current_frame->r[j].uz - temp_frame->r[j].uz;
524    
525     dx /= (double)(n_interpolate + 1);
526     dy /= (double)(n_interpolate + 1);
527     dz /= (double)(n_interpolate + 1);
528    
529     out_coords[j].hasVector = current_frame->r[j].hasVector;
530     out_coords[j].ux = temp_frame->r[j].ux + dx * (i+1);
531     out_coords[j].uy = temp_frame->r[j].uy + dy * (i+1);
532     out_coords[j].uz = temp_frame->r[j].uz + dz * (i+1);
533     }
534    
535     if (current_frame->r[j].hasCharge) {
536     dx = current_frame->r[j].charge - temp_frame->r[j].charge;
537    
538     dx /= (double)(n_interpolate + 1);
539    
540     out_coords[j].hasCharge = current_frame->r[j].hasCharge;
541     out_coords[j].charge = temp_frame->r[j].charge + dx * (i+1);
542     }
543    
544 mmeineke 515 }
545    
546     pov_write(out_file, out_coords, n_atoms, draw_hydrogens, draw_bonds,
547 gezelter 864 draw_atoms, draw_vectors);
548 mmeineke 515 free(out_coords);
549     (void)fclose(out_file);
550 mmeineke 60 }
551 mmeineke 515 }
552    
553     /* open the new output file */
554    
555     sprintf(out_name, out_format, currentCount );
556     out_file = fopen(out_name, "w");
557     if(out_file == NULL){
558     printf("error opening output file: %s\n", out_name);
559     exit(8);
560     }
561     (void)fprintf(out_file,
562     "// The following script was automatically generated by:\n"
563     "// xyz2pov Copyright 2001 by MATTHEW A. MEINEKE\n"
564     "\n"
565     "#include \"pov_header.pov\"\n"
566     "\n");
567    
568     if( draw_box ){
569 mmeineke 60
570 mmeineke 515 fprintf( out_file,
571 gezelter 630 "makePeriodicBox( %lf, %lf, %lf, %lf, %lf, %lf, %lf, %lf, %lf )\n"
572 mmeineke 515 "\n",
573 gezelter 625 current_frame->Hmat[0][0],
574 mmeineke 629 current_frame->Hmat[2][0],
575 gezelter 625 current_frame->Hmat[1][0],
576     current_frame->Hmat[0][1],
577 mmeineke 629 current_frame->Hmat[2][1],
578 gezelter 625 current_frame->Hmat[1][1],
579     current_frame->Hmat[0][2],
580 mmeineke 629 current_frame->Hmat[2][2],
581     current_frame->Hmat[1][2] );
582 mmeineke 60 }
583 mmeineke 515
584    
585    
586     out_coords =
587     (struct coords *)calloc(n_atoms, sizeof(struct coords));
588    
589     for(i = 0; i < n_atoms; i++){
590     strcpy(out_coords[i].name, current_frame->r[i].name);
591     out_coords[i].x = current_frame->r[i].x;
592     out_coords[i].y = current_frame->r[i].y;
593     out_coords[i].z = current_frame->r[i].z;
594 gezelter 864
595     if (current_frame->r[i].hasVector) {
596     out_coords[i].hasVector = current_frame->r[i].hasVector;
597     out_coords[i].ux = current_frame->r[i].ux;
598     out_coords[i].uy = current_frame->r[i].uy;
599     out_coords[i].uz = current_frame->r[i].uz;
600     }
601    
602     if (current_frame->r[i].hasCharge) {
603     out_coords[i].hasCharge = current_frame->r[i].hasCharge;
604     out_coords[i].charge = current_frame->r[i].charge;
605     }
606 mmeineke 515 }
607     pov_write(out_file, out_coords, n_atoms, draw_hydrogens, draw_bonds,
608 gezelter 864 draw_atoms, draw_vectors);
609 mmeineke 515 free(out_coords);
610    
611     (void)fclose(out_file);
612 mmeineke 60 }
613    
614     /*free up memory */
615    
616     temp_frame = current_frame->next;
617     current_frame->next = NULL;
618    
619     if(temp_frame != NULL){
620    
621     free(temp_frame->r);
622     free(temp_frame);
623     }
624    
625     /* make a new frame */
626    
627     temp_frame = (struct linked_xyz *)malloc(sizeof(struct linked_xyz));
628     temp_frame->next = current_frame;
629     current_frame = temp_frame;
630    
631     eof_test = fgets(read_buffer, sizeof(read_buffer), in_file);
632    
633 mmeineke 515 if( haveEnd ){
634     if( currentCount >= (endFrame + n_interpolate * span) ) done = 1;
635     }
636 mmeineke 60 }
637    
638     (void)fclose(in_file);
639    
640    
641     if(generate_header){
642    
643     dx = big_x - small_x;
644     dy = big_y - small_y;
645     dz = big_z - small_z;
646    
647     rsqr = dx * dx + dy * dy + dz * dz;
648     diagonal = sqrt(rsqr);
649     diagonal *= 0.5;
650    
651 mmeineke 515 // calculate the center
652    
653     dx = big_x + small_x;
654     dy = big_y + small_y;
655     dz = big_z + small_z;
656    
657 mmeineke 60 dx /= 2.0;
658     dy /= 2.0;
659     dz /= 2.0;
660 mmeineke 515
661 mmeineke 60
662     /*note the y and z axis is exchanged for the different coordinate
663     system in pov-ray*/
664    
665    
666     out_file = fopen("pov_header.pov", "w");
667    
668     fprintf(out_file,
669     "// The following script was automatically generated by:\n"
670     "// xyz2pov Copyright 2001 by MATTHEW A. MEINEKE\n"
671     "\n"
672     "\n"
673     "background { rgb <1.0, 1.0, 1.0> }\n"
674     "\n"
675     "\n"
676     );
677    
678     fprintf(out_file,
679     "//******************************************************\n"
680     "// Declare the resolution, camera, and light sources.\n"
681     "//******************************************************\n"
682     "\n"
683     "// NOTE: if you plan to render at a different resoltion,\n"
684     "// be sure to update the following two lines to maintain\n"
685     "// the correct aspect ratio.\n"
686     "\n"
687     "#declare Width = 640.0;\n"
688     "#declare Height = 480.0;\n"
689     "#declare Ratio = Width / Height;\n"
690     "\n"
691     "#declare ATOM_SPHERE_FACTOR = 0.2;\n"
692     "#declare BOND_RADIUS = 0.1;\n"
693 gezelter 864 "#declare VECTOR_SCALE = 1.0;\n"
694     "#declare STICK_RADIUS = 0.5 * BOND_RADIUS;\n"
695     "#declare CONE_RADIUS = 2.0 * STICK_RADIUS;\n"
696     "#declare CONE_FRACTION = 0.15;\n"
697 mmeineke 60 "\n"
698 mmeineke 515 "// declare camera, light, and system variables\n"
699     "\n"
700     "#declare sysCenterX = %lf;\n"
701     "#declare sysCenterY = %lf;\n"
702     "#declare sysCenterZ = %lf;\n"
703     "\n"
704     "#declare zoom = %lf;\n"
705     "\n",
706     dx, dz, dy,
707     diagonal );
708    
709     fprintf(out_file,
710     "#declare cameraLookX = sysCenterX;\n"
711     "#declare cameraLookY = sysCenterY;\n"
712     "#declare cameraLookZ = sysCenterZ;\n"
713     "\n"
714 gezelter 868 "#declare rotatePointX = cameraLookX;\n"
715     "#declare rotatePointY = cameraLookY;\n"
716     "#declare rotatePointZ = cameraLookZ;\n"
717     "\n"
718 mmeineke 515 "#declare cameraX = cameraLookX;\n"
719     "#declare cameraY = cameraLookY;\n"
720     "#declare cameraZ = cameraLookZ - zoom;\n"
721     "\n"
722     "#declare lightAx = cameraX;\n"
723     "#declare lightAy = cameraY;\n"
724     "#declare lightAz = cameraZ;\n"
725     "\n"
726     "#declare lightBx = cameraX - zoom;\n"
727     "#declare lightBy = cameraY + zoom;\n"
728     "#declare lightBz = cameraZ;\n"
729     "\n"
730     "#declare boxCenterX = cameraLookX;\n"
731     "#declare boxCenterY = cameraLookY;\n"
732     "#declare boxCenterZ = cameraLookZ;\n"
733     "\n"
734     "// declare the cameras and the light sources\n"
735     "\n"
736 mmeineke 60 "camera{\n"
737 mmeineke 515 " location < cameraX, cameraY, cameraZ>\n"
738 mmeineke 60 " right < Ratio , 0, 0>\n"
739 mmeineke 515 " look_at < cameraLookX, cameraLookY, cameraLookZ >\n"
740 mmeineke 60 "}\n"
741 mmeineke 515 "\n"
742 mmeineke 60 "light_source{\n"
743 mmeineke 515 " < lightAx, lightAy, lightAz >\n"
744     " rgb < 1.0, 1.0, 1.0 > }\n"
745     "\n"
746 mmeineke 60 "light_source{\n"
747 mmeineke 515 " < lightBx, lightBy, lightBz >\n"
748 mmeineke 60 " rgb < 1.0, 1.0, 1.0 > }\n"
749     "\n"
750 mmeineke 515 "\n"
751 mmeineke 60 "//************************************************************\n"
752     "// Set whether or not to rotate the system.\n"
753     "//\n"
754     "// To Rotate, set ROTATE to 1.0 (true),\n"
755     "// Then set the Euler Angles PHI, THETA, and PSI in degrees.\n"
756     "//************************************************************\n"
757     "\n"
758     "#declare ROTATE = 0.0;\n"
759     "#declare PHI = 0.0;\n"
760     "#declare THETA = 0.0;\n"
761     "#declare PSI = 0.0;\n"
762     "\n"
763     "#if(ROTATE)\n"
764     " #declare phi_r = radians(PHI);\n"
765     " #declare theta_r = radians(THETA);\n"
766     " #declare psi_r = radians(PSI);\n"
767     "\n"
768     " #declare A11 = cos(phi_r) * cos(psi_r) - sin(phi_r) * cos(theta_r) * sin(psi_r);\n"
769     " #declare A12 = sin(phi_r) * cos(psi_r) + cos(phi_r) * cos(theta_r) * sin(psi_r);\n"
770     " #declare A13 = sin(theta_r) * sin(psi_r);\n"
771     "\n"
772     " #declare A21 = -cos(phi_r) * sin(psi_r) - sin(phi_r) * cos(theta_r) * cos(psi_r);\n"
773     " #declare A22 = -sin(phi_r) * sin(psi_r) + cos(phi_r) * cos(theta_r) * cos(psi_r);\n"
774     " #declare A23 = sin(theta_r) * cos(psi_r);\n"
775     "\n"
776     " #declare A31 = sin(phi_r) * sin(theta_r);\n"
777     " #declare A32 = -cos(phi_r) * sin(theta_r);\n"
778     " #declare A33 = cos(theta_r);\n"
779     "\n"
780     "#end\n"
781 mmeineke 515 "\n"
782     "\n"
783     "//************************************************************\n"
784     "// declare the periodic box macro\n"
785     "//************************************************************\n"
786     "\n"
787 gezelter 625 "#macro makePeriodicBox( bx1, by1, bz1, bx2, by2, bz2, bx3, by3, bz3 )\n"
788 mmeineke 515 "\n"
789 mmeineke 629 " #local bcx = (bx1 + bx2 + bx3) / 2.0;\n"
790     " #local bcy = (by1 + by2 + by3) / 2.0;\n"
791     " #local bcz = (bz1 + bz2 + bz3) / 2.0;\n"
792 mmeineke 515 "\n"
793 mmeineke 629 " #local pAx = boxCenterX - bcx;\n"
794     " #local pAy = boxCenterY - bcy;\n"
795     " #local pAz = boxCenterZ - bcz;\n"
796     " #local pBx = boxCenterX + bx1 - bcx;\n"
797     " #local pBy = boxCenterY + by1 - bcy;\n"
798     " #local pBz = boxCenterZ + bz1 - bcz;\n"
799     " #local pCx = boxCenterX + bx2 - bcx;\n"
800     " #local pCy = boxCenterY + by2 - bcy;\n"
801     " #local pCz = boxCenterZ + bz2 - bcz;\n"
802     " #local pDx = boxCenterX + bx3 - bcx;\n"
803     " #local pDy = boxCenterY + by3 - bcy;\n"
804     " #local pDz = boxCenterZ + bz3 - bcz;\n"
805     " #local pEx = boxCenterX + bx1 + bx2 - bcx;\n"
806     " #local pEy = boxCenterY + by1 + by2 - bcy;\n"
807     " #local pEz = boxCenterZ + bz1 + bz2 - bcz;\n"
808     " #local pFx = boxCenterX + bx1 + bx3 - bcx;\n"
809     " #local pFy = boxCenterY + by1 + by3 - bcy;\n"
810     " #local pFz = boxCenterZ + bz1 + bz3 - bcz;\n"
811     " #local pGx = boxCenterX + bx2 + bx3 - bcx;\n"
812     " #local pGy = boxCenterY + by2 + by3 - bcy;\n"
813     " #local pGz = boxCenterZ + bz2 + bz3 - bcz;\n"
814     " #local pHx = boxCenterX + bx1 + bx2 + bx3 - bcx;\n"
815     " #local pHy = boxCenterY + by1 + by2 + by3 - bcy;\n"
816     " #local pHz = boxCenterZ + bz1 + bz2 + bz3 - bcz;\n"
817 gezelter 868 "\n"
818     " #if(ROTATE)\n"
819     " #local pAx_new = rotatePointX + A11 * (pAx-rotatePointX) + A12 * (pAy-rotatePointY) + A13 * (pAz-rotatePointZ);\n"
820     " #local pAy_new = rotatePointY + A21 * (pAx-rotatePointX) + A22 * (pAy-rotatePointY) + A23 * (pAz-rotatePointZ);\n"
821     " #local pAz_new = rotatePointZ + A31 * (pAx-rotatePointX) + A32 * (pAy-rotatePointY) + A33 * (pAz-rotatePointZ);\n"
822 mmeineke 629 "\n"
823 gezelter 868 " #local pBx_new = rotatePointX + A11 * (pBx-rotatePointX) + A12 * (pBy-rotatePointY) + A13 * (pBz-rotatePointZ);\n"
824     " #local pBy_new = rotatePointY + A21 * (pBx-rotatePointX) + A22 * (pBy-rotatePointY) + A23 * (pBz-rotatePointZ);\n"
825     " #local pBz_new = rotatePointZ + A31 * (pBx-rotatePointX) + A32 * (pBy-rotatePointY) + A33 * (pBz-rotatePointZ);\n"
826     "\n"
827     " #local pCx_new = rotatePointX + A11 * (pCx-rotatePointX) + A12 * (pCy-rotatePointY) + A13 * (pCz-rotatePointZ);\n"
828     " #local pCy_new = rotatePointY + A21 * (pCx-rotatePointX) + A22 * (pCy-rotatePointY) + A23 * (pCz-rotatePointZ);\n"
829     " #local pCz_new = rotatePointZ + A31 * (pCx-rotatePointX) + A32 * (pCy-rotatePointY) + A33 * (pCz-rotatePointZ);\n"
830     "\n"
831     " #local pDx_new = rotatePointX + A11 * (pDx-rotatePointX) + A12 * (pDy-rotatePointY) + A13 * (pDz-rotatePointZ);\n"
832     " #local pDy_new = rotatePointY + A21 * (pDx-rotatePointX) + A22 * (pDy-rotatePointY) + A23 * (pDz-rotatePointZ);\n"
833     " #local pDz_new = rotatePointZ + A31 * (pDx-rotatePointX) + A32 * (pDy-rotatePointY) + A33 * (pDz-rotatePointZ);\n"
834     "\n"
835     " #local pEx_new = rotatePointX + A11 * (pEx-rotatePointX) + A12 * (pEy-rotatePointY) + A13 * (pEz-rotatePointZ);\n"
836     " #local pEy_new = rotatePointY + A21 * (pEx-rotatePointX) + A22 * (pEy-rotatePointY) + A23 * (pEz-rotatePointZ);\n"
837     " #local pEz_new = rotatePointZ + A31 * (pEx-rotatePointX) + A32 * (pEy-rotatePointY) + A33 * (pEz-rotatePointZ);\n"
838     "\n"
839     " #local pFx_new = rotatePointX + A11 * (pFx-rotatePointX) + A12 * (pFy-rotatePointY) + A13 * (pFz-rotatePointZ);\n"
840     " #local pFy_new = rotatePointY + A21 * (pFx-rotatePointX) + A22 * (pFy-rotatePointY) + A23 * (pFz-rotatePointZ);\n"
841     " #local pFz_new = rotatePointZ + A31 * (pFx-rotatePointX) + A32 * (pFy-rotatePointY) + A33 * (pFz-rotatePointZ);\n"
842     "\n"
843     " #local pGx_new = rotatePointX + A11 * (pGx-rotatePointX) + A12 * (pGy-rotatePointY) + A13 * (pGz-rotatePointZ);\n"
844     " #local pGy_new = rotatePointY + A21 * (pGx-rotatePointX) + A22 * (pGy-rotatePointY) + A23 * (pGz-rotatePointZ);\n"
845     " #local pGz_new = rotatePointZ + A31 * (pGx-rotatePointX) + A32 * (pGy-rotatePointY) + A33 * (pGz-rotatePointZ);\n"
846     "\n"
847     " #local pHx_new = rotatePointX + A11 * (pHx-rotatePointX) + A12 * (pHy-rotatePointY) + A13 * (pHz-rotatePointZ);\n"
848     " #local pHy_new = rotatePointY + A21 * (pHx-rotatePointX) + A22 * (pHy-rotatePointY) + A23 * (pHz-rotatePointZ);\n"
849     " #local pHz_new = rotatePointZ + A31 * (pHx-rotatePointX) + A32 * (pHy-rotatePointY) + A33 * (pHz-rotatePointZ);\n"
850     "\n"
851     " #else\n"
852     " #local pAx_new = pAx;"
853     " #local pAy_new = pAy;"
854     " #local pAz_new = pAz;"
855     "\n"
856     " #local pBx_new = pBx;"
857     " #local pBy_new = pBy;"
858     " #local pBz_new = pBz;"
859     "\n"
860     " #local pCx_new = pCx;"
861     " #local pCy_new = pCy;"
862     " #local pCz_new = pCz;"
863     "\n"
864     " #local pDx_new = pDx;"
865     " #local pDy_new = pDy;"
866     " #local pDz_new = pDz;"
867     "\n"
868     " #local pEx_new = pEx;"
869     " #local pEy_new = pEy;"
870     " #local pEz_new = pEz;"
871     "\n"
872     " #local pFx_new = pFx;"
873     " #local pFy_new = pFy;"
874     " #local pFz_new = pFz;"
875     "\n"
876     " #local pGx_new = pGx;"
877     " #local pGy_new = pGy;"
878     " #local pGz_new = pGz;"
879     "\n"
880     " #local pHx_new = pHx;"
881     " #local pHy_new = pHy;"
882     " #local pHz_new = pHz;"
883     "\n"
884     " #end\n"
885     " #local pAx = pAx_new;"
886     " #local pAy = pAy_new;"
887     " #local pAz = pAz_new;"
888     "\n"
889     " #local pBx = pBx_new;"
890     " #local pBy = pBy_new;"
891     " #local pBz = pBz_new;"
892     "\n"
893     " #local pCx = pCx_new;"
894     " #local pCy = pCy_new;"
895     " #local pCz = pCz_new;"
896     "\n"
897     " #local pDx = pDx_new;"
898     " #local pDy = pDy_new;"
899     " #local pDz = pDz_new;"
900     "\n"
901     " #local pEx = pEx_new;"
902     " #local pEy = pEy_new;"
903     " #local pEz = pEz_new;"
904     "\n"
905     " #local pFx = pFx_new;"
906     " #local pFy = pFy_new;"
907     " #local pFz = pFz_new;"
908     "\n"
909     " #local pGx = pGx_new;"
910     " #local pGy = pGy_new;"
911     " #local pGz = pGz_new;"
912     "\n"
913     " #local pHx = pHx_new;"
914     " #local pHy = pHy_new;"
915     " #local pHz = pHz_new;"
916     "\n"
917 mmeineke 515 " #local colorR = 0.90;\n"
918     " #local colorG = 0.91;\n"
919     " #local colorB = 0.98;\n"
920     "\n"
921     " #local pipeWidth = 0.4;\n"
922     "\n"
923     " // 1\n"
924     " cylinder{\n"
925 gezelter 625 " < pAx, pAy, pAz >,\n"
926     " < pBx, pBy, pBz >,\n"
927 mmeineke 515 " pipeWidth\n"
928     " texture{\n"
929     " pigment{ rgb < colorR, colorG, colorB > }\n"
930     " finish{\n"
931     " ambient .2\n"
932     " diffuse .6\n"
933     " specular 1\n"
934     " roughness .001\n"
935     " metallic\n"
936     " }\n"
937     " }\n"
938     " }\n"
939     "\n"
940     " // 2\n"
941     " cylinder{\n"
942 gezelter 625 " < pAx, pAy, pAz >,\n"
943     " < pCx, pCy, pCz >,\n"
944 mmeineke 515 " pipeWidth\n"
945     " texture{\n"
946     " pigment{ rgb < colorR, colorG, colorB > }\n"
947     " finish{\n"
948     " ambient .2\n"
949     " diffuse .6\n"
950     " specular 1\n"
951     " roughness .001\n"
952     " metallic\n"
953     " }\n"
954     " }\n"
955     " }\n"
956     "\n"
957     " // 3\n"
958     " cylinder{\n"
959 gezelter 625 " < pAx, pAy, pAz >,\n"
960     " < pDx, pDy, pDz >,\n"
961 mmeineke 515 " pipeWidth\n"
962     " texture{\n"
963     " pigment{ rgb < colorR, colorG, colorB > }\n"
964     " finish{\n"
965     " ambient .2\n"
966     " diffuse .6\n"
967     " specular 1\n"
968     " roughness .001\n"
969     " metallic\n"
970     " }\n"
971     " }\n"
972     " }\n"
973     "\n"
974     " // 4\n"
975     " cylinder{\n"
976 gezelter 625 " < pBx, pBy, pBz >,\n"
977     " < pEx, pEy, pEz >,\n"
978 mmeineke 515 " pipeWidth\n"
979     " texture{\n"
980     " pigment{ rgb < colorR, colorG, colorB > }\n"
981     " finish{\n"
982     " ambient .2\n"
983     " diffuse .6\n"
984     " specular 1\n"
985     " roughness .001\n"
986     " metallic\n"
987     " }\n"
988     " }\n"
989     " }\n"
990     "\n"
991     " // 5\n"
992     " cylinder{\n"
993 gezelter 625 " < pCx, pCy, pCz >,\n"
994     " < pEx, pEy, pEz >,\n"
995 mmeineke 515 " pipeWidth\n"
996     " texture{\n"
997     " pigment{ rgb < colorR, colorG, colorB > }\n"
998     " finish{\n"
999     " ambient .2\n"
1000     " diffuse .6\n"
1001     " specular 1\n"
1002     " roughness .001\n"
1003     " metallic\n"
1004     " }\n"
1005     " }\n"
1006     " }\n"
1007     "\n"
1008     " // 6\n"
1009     " cylinder{\n"
1010 gezelter 625 " < pBx, pBy, pBz >,\n"
1011     " < pFx, pFy, pFz >,\n"
1012 mmeineke 515 " pipeWidth\n"
1013     " texture{\n"
1014     " pigment{ rgb < colorR, colorG, colorB > }\n"
1015     " finish{\n"
1016     " ambient .2\n"
1017     " diffuse .6\n"
1018     " specular 1\n"
1019     " roughness .001\n"
1020     " metallic\n"
1021     " }\n"
1022     " }\n"
1023     " }\n"
1024     "\n"
1025     " // 7\n"
1026     " cylinder{\n"
1027 gezelter 625 " < pCx, pCy, pCz >,\n"
1028     " < pGx, pGy, pGz >,\n"
1029 mmeineke 515 " pipeWidth\n"
1030     " texture{\n"
1031     " pigment{ rgb < colorR, colorG, colorB > }\n"
1032     " finish{\n"
1033     " ambient .2\n"
1034     " diffuse .6\n"
1035     " specular 1\n"
1036     " roughness .001\n"
1037     " metallic\n"
1038     " }\n"
1039     " }\n"
1040     " }\n"
1041     "\n"
1042     " // 8\n"
1043     " cylinder{\n"
1044 gezelter 625 " < pDx, pDy, pDz >,\n"
1045     " < pGx, pGy, pGz >,\n"
1046 mmeineke 515 " pipeWidth\n"
1047     " texture{\n"
1048     " pigment{ rgb < colorR, colorG, colorB > }\n"
1049     " finish{\n"
1050     " ambient .2\n"
1051     " diffuse .6\n"
1052     " specular 1\n"
1053     " roughness .001\n"
1054     " metallic\n"
1055     " }\n"
1056     " }\n"
1057     " }\n"
1058     "\n"
1059     " // 9\n"
1060     " cylinder{\n"
1061 gezelter 625 " < pDx, pDy, pDz >,\n"
1062     " < pFx, pFy, pFz >,\n"
1063 mmeineke 515 " pipeWidth\n"
1064     " texture{\n"
1065     " pigment{ rgb < colorR, colorG, colorB > }\n"
1066     " finish{\n"
1067     " ambient .2\n"
1068     " diffuse .6\n"
1069     " specular 1\n"
1070     " roughness .001\n"
1071     " metallic\n"
1072     " }\n"
1073     " }\n"
1074     " }\n"
1075     "\n"
1076     " // 10\n"
1077     " cylinder{\n"
1078 gezelter 625 " < pEx, pEy, pEz >,\n"
1079     " < pHx, pHy, pHz >,\n"
1080 mmeineke 515 " pipeWidth\n"
1081     " texture{\n"
1082     " pigment{ rgb < colorR, colorG, colorB > }\n"
1083     " finish{\n"
1084     " ambient .2\n"
1085     " diffuse .6\n"
1086     " specular 1\n"
1087     " roughness .001\n"
1088     " metallic\n"
1089     " }\n"
1090     " }\n"
1091     " }\n"
1092     "\n"
1093     " // 11\n"
1094     " cylinder{\n"
1095 gezelter 625 " < pFx, pFy, pFz >,\n"
1096     " < pHx, pHy, pHz >,\n"
1097 mmeineke 515 " pipeWidth\n"
1098     " texture{\n"
1099     " pigment{ rgb < colorR, colorG, colorB > }\n"
1100     " finish{\n"
1101     " ambient .2\n"
1102     " diffuse .6\n"
1103     " specular 1\n"
1104     " roughness .001\n"
1105     " metallic\n"
1106     " }\n"
1107     " }\n"
1108     " }\n"
1109     "\n"
1110     " // 12\n"
1111     " cylinder{\n"
1112 gezelter 625 " < pGx, pGy, pGz >,\n"
1113     " < pHx, pHy, pHz >,\n"
1114 mmeineke 515 " pipeWidth\n"
1115     " texture{\n"
1116     " pigment{ rgb < colorR, colorG, colorB > }\n"
1117     " finish{\n"
1118     " ambient .2\n"
1119     " diffuse .6\n"
1120     " specular 1\n"
1121     " roughness .001\n"
1122     " metallic\n"
1123     " }\n"
1124     " }\n"
1125     " }\n"
1126     "\n"
1127     "#end\n"
1128     "\n"
1129 mmeineke 60 "\n");
1130 mmeineke 515
1131    
1132    
1133    
1134 mmeineke 60 make_header_macros(out_file);
1135    
1136     fclose(out_file);
1137     }
1138    
1139     return 0;
1140    
1141     }
1142    
1143    
1144    
1145     /***************************************************************************
1146     * prints out the usage for the command line arguments, then exits.
1147     ***************************************************************************/
1148    
1149     void usage(){
1150     (void)fprintf(stderr,
1151     "The proper usage is: %s [options] <xyz_file>\n"
1152     "\n"
1153     "Options:\n"
1154     "\n"
1155     " -o <prefix> the output file prefix\n"
1156     " -H generate a pov-ray header file\n"
1157     " -i <#> number of frames to interpolate\n"
1158     " -h draw hydrogens\n"
1159     " -b draw bonds\n"
1160     " -a draw atoms\n"
1161 gezelter 864 " -v draw vectors\n"
1162 mmeineke 515 " -p draw periodic box\n"
1163 mmeineke 508 " -r regenerate bond\n"
1164 mmeineke 515 " -f <#> render frame <#> only\n"
1165     " -s <#> begin at frame <#> (inclusive)\n"
1166     " -e <#> end at frame <#> (inclusive)\n"
1167 mmeineke 60 "\n",
1168     program_name);
1169     exit(8);
1170     }
1171 gezelter 864
1172     int count_tokens(line, delimiters)
1173     /* PURPOSE: RETURN A COUNT OF THE NUMBER OF TOKENS ON THE LINE. */
1174     char *line; /* LINE CONTAINING TOKENS. */
1175     char *delimiters; /* POSSIBLE TOKEN DELIMITERS TO USE. */
1176     {
1177     char *working_line; /* WORKING COPY OF LINE. */
1178     int ntokens; /* NUMBER OF TOKENS FOUND IN LINE. */
1179     char *strtok_ptr; /* POINTER FOR STRTOK. */
1180    
1181     strtok_ptr= working_line= strdup(line);
1182    
1183     ntokens=0;
1184     while (strtok(strtok_ptr,delimiters)!=NULL)
1185     {
1186     ntokens++;
1187     strtok_ptr=NULL;
1188     }
1189    
1190     free(working_line);
1191     return(ntokens);
1192     }