ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/group/trunk/madProps/madProps.c
Revision: 82
Committed: Fri Aug 16 04:05:45 2002 UTC (21 years, 11 months ago) by mmeineke
Content type: text/plain
File size: 9569 byte(s)
Log Message:
added the capability to specify a start and end frame for cosCorr and GofR

File Contents

# User Rev Content
1 mmeineke 38 #include <stdio.h>
2     #include <stdlib.h>
3     #include <string.h>
4     #include <math.h>
5    
6 mmeineke 43
7     #include "madProps.h"
8 mmeineke 38 #include "frameCount.h"
9    
10     char *program_name; /*the name of the program */
11    
12     void usage(void);
13    
14     int main(argc, argv)
15     int argc;
16     char *argv[];
17     {
18    
19    
20     struct xyz_frame *dumpArray;
21     int nFrames; // the nu8mvber of frames
22     int lineNum = 0; // keeps track of the line number
23     int n_atoms; // the number of atoms
24     int i,j; // loop counters
25 mmeineke 45 int isFirst;
26 mmeineke 38
27     char read_buffer[2000]; /*the line buffer for reading */
28     char *foo; /*the pointer to the current string token */
29     FILE *in_file; /* the input file */
30     char *in_name = NULL; /*the name of the input file */
31 mmeineke 43 char *out_prefix; // the output prefix
32     char current_flag; // used in parseing the flags
33 mmeineke 38
34    
35 mmeineke 43
36     int done = 0; // multipurpose boolean
37     int have_prefix = 0; // prefix boolean
38     int calcRMSD = 0;
39    
40     int calcGofR = 0;
41     char gofR1[30];
42     char gofR2[30];
43    
44     int calcMuCorr = 0;
45     char muCorr[30];
46    
47     int calcCosCorr = 0;
48     char cosCorr1[30];
49     char cosCorr2[30];
50    
51 mmeineke 82 int startFrame = 0;
52     int haveStartFrame = 0;
53     int endFrame = 0;
54     int haveEndFrame = 0;
55    
56 mmeineke 38 program_name = argv[0]; /*save the program name in case we need it*/
57    
58 mmeineke 43 for( i = 1; i < argc; i++){
59 mmeineke 38
60     if(argv[i][0] =='-'){
61 mmeineke 43
62     // parse the option
63 mmeineke 38
64 mmeineke 43 if(argv[i][1] == '-' ){
65    
66     // parse long word options
67    
68     if( !strcmp( argv[i], "--GofR" ) ){
69     calcGofR = 1;
70     i++;
71     strcpy( gofR1, argv[i] );
72     i++;
73     strcpy( gofR2, argv[i] );
74     }
75    
76     else if( !strcmp( argv[i], "--MuCorr") ){
77     calcMuCorr = 1;
78     i++;
79     strcpy( muCorr, argv[i] );
80     }
81    
82 mmeineke 82 else if( !strcmp( argv[i], "--startFrame" ) ){
83     haveStartFrame = 1;
84 mmeineke 43 i++;
85 mmeineke 82 startFrame = atoi(argv[i]);
86     }
87    
88     else if( !strcmp( argv[i], "--endFrame" ) ){
89     haveEndFrame = 1;
90 mmeineke 43 i++;
91 mmeineke 82 endFrame = atoi(argv[i]);
92 mmeineke 43 }
93    
94     else{
95     fprintf( stderr,
96     "Invalid option \"%s\"\n", argv[i] );
97     usage();
98     }
99     }
100 mmeineke 38
101 mmeineke 43 else{
102 mmeineke 38
103 mmeineke 43 // parse single character options
104    
105     done =0;
106     j = 1;
107     current_flag = argv[i][j];
108     while( (current_flag != '\0') && (!done) ){
109    
110     switch(current_flag){
111 mmeineke 38
112 mmeineke 43 case 'o':
113     // -o <prefix> => the output prefix.
114 mmeineke 38
115 mmeineke 43 i++;
116     out_prefix = argv[i];
117     have_prefix = 1;
118     done = 1;
119     break;
120    
121     case 'h':
122     // -h => give the usage
123    
124     usage();
125     break;
126    
127     case 'r':
128     // calculates the rmsd
129    
130     calcRMSD = 1;
131     break;
132    
133     case 'g':
134     // calculate all to all g(r)
135    
136     calcGofR = 1;
137     strcpy( gofR1, "all" );
138     strcpy( gofR2, "all" );
139     break;
140    
141     default:
142    
143     fprintf( stderr, "about to print bad option\n" );
144    
145     (void)fprintf(stderr, "Bad option \"-%s\"\n", current_flag);
146     usage();
147     }
148     j++;
149     current_flag = argv[i][j];
150     }
151     }
152     }
153    
154     else{
155    
156     if( in_name != NULL ){
157     fprintf( stderr,
158     "Error at \"%s\", program does not currently support\n"
159     "more than one input file.\n"
160     "\n",
161     argv[i]);
162 mmeineke 38 usage();
163     }
164 mmeineke 43
165     in_name = argv[i];
166 mmeineke 38 }
167     }
168    
169     if(in_name == NULL){
170     usage();
171     }
172 mmeineke 43
173     if( !have_prefix ) out_prefix = in_name;
174 mmeineke 38
175     printf( "Counting number of frames..." );
176     fflush( stdout );
177    
178     nFrames = frameCount( in_name );
179 mmeineke 82 if( !haveEndFrame ) endFrame = nFrames;
180 mmeineke 38
181     printf( "done.\n"
182     "nframes = %d\n"
183     "\n",
184     nFrames );
185     fflush( stdout );
186    
187     in_file = fopen(in_name, "r");
188     if(in_file == NULL){
189     printf("Cannot open file: %s\n", in_name);
190     exit(8);
191     }
192    
193 mmeineke 45 // create and initialize the array of frames
194 mmeineke 38
195     dumpArray = (struct xyz_frame*)calloc( nFrames,
196     sizeof( struct xyz_frame ) );
197 mmeineke 45 for( i=0; i<nFrames; i++ ){
198     dumpArray[i].nAtoms = 0;
199     dumpArray[i].time = 0.0;
200     dumpArray[i].boxX = 0.0;
201     dumpArray[i].boxY = 0.0;
202     dumpArray[i].boxZ = 0.0;
203     dumpArray[i].r = NULL;
204     dumpArray[i].v = NULL;
205     dumpArray[i].names = NULL;
206     }
207 mmeineke 38
208     // read the frames
209    
210     printf( "Reading the frames into the coordinate arrays..." );
211     fflush( stdout );
212    
213 mmeineke 45 isFirst = 1;
214 mmeineke 38 for(j =0; j<nFrames; j++ ){
215    
216     // read the number of atoms
217    
218     fgets(read_buffer, sizeof(read_buffer), in_file);
219     lineNum++;
220    
221     n_atoms = atoi( read_buffer );
222     dumpArray[j].nAtoms = n_atoms;
223    
224     dumpArray[j].r =
225     (struct coords *)calloc(n_atoms, sizeof(struct coords));
226 mmeineke 45
227     if( isFirst ) {
228     dumpArray[0].names =
229     (atomID *)calloc( n_atoms, sizeof(atomID) );
230     isFirst = 0;
231     }
232 mmeineke 38
233 mmeineke 43 if( calcMuCorr || calcCosCorr ){
234     dumpArray[j].v =
235     (struct vect *)calloc(n_atoms, sizeof(struct vect));
236     }
237    
238 mmeineke 38 //read the time and the box sizes
239    
240     fgets(read_buffer, sizeof(read_buffer), in_file);
241     lineNum++;
242    
243     foo = strtok(read_buffer, " ,;\t");
244     if(foo == NULL){
245     printf("error in reading time at line %d\n", lineNum);
246     exit(8);
247     }
248    
249     dumpArray[j].time = atof( foo );
250    
251     foo = strtok(NULL, " ,;\t");
252     if(foo == NULL){
253     printf("error in reading boxX at line %d\n", lineNum);
254     exit(8);
255     }
256    
257     dumpArray[j].boxX = atof( foo );
258    
259     foo = strtok(NULL, " ,;\t");
260     if(foo == NULL){
261     printf("error in reading boxY at line %d\n", lineNum);
262     exit(8);
263     }
264    
265     dumpArray[j].boxY = atof( foo );
266    
267     foo = strtok(NULL, " ,;\t");
268     if(foo == NULL){
269     printf("error in reading boxZ at line %d\n", lineNum);
270     exit(8);
271     }
272    
273     dumpArray[j].boxZ = atof( foo );
274    
275    
276     for( i=0; i < n_atoms; i++){
277    
278     fgets(read_buffer, sizeof(read_buffer), in_file);
279     lineNum++;
280    
281     // get the name and positions
282    
283     foo = strtok(read_buffer, " ,;\t");
284     if(foo == NULL){
285     printf("error in reading atom name at line %d\n", lineNum );
286     exit(8);
287     }
288    
289 mmeineke 45 strcpy(dumpArray[0].names[i], foo); /*copy the atom name */
290 mmeineke 38
291     foo = strtok(NULL, " ,;\t");
292     if(foo == NULL){
293     printf("error in reading position x at line %d\n", lineNum);
294     exit(8);
295     }
296    
297     dumpArray[j].r[i].x = atof( foo );
298    
299     foo = strtok(NULL, " ,;\t");
300     if(foo == NULL){
301     printf("error in reading position y at line %d\n", lineNum);
302     exit(8);
303     }
304    
305     dumpArray[j].r[i].y = atof( foo );
306    
307     foo = strtok(NULL, " ,;\t");
308     if(foo == NULL){
309     printf("error in reading position z at line %d\n", lineNum);
310     exit(8);
311     }
312    
313     dumpArray[j].r[i].z = atof( foo );
314 mmeineke 43
315     if( calcCosCorr || calcMuCorr ){
316    
317     foo = strtok(NULL, " ,;\t");
318     if(foo == NULL){
319    
320     dumpArray[j].v[i].x = 0.0;
321     dumpArray[j].v[i].y = 0.0;
322     dumpArray[j].v[i].z = 0.0;
323     }
324     else{
325    
326     dumpArray[j].v[i].x = atof( foo );
327    
328     foo = strtok(NULL, " ,;\t");
329     if(foo == NULL){
330     printf("error in reading vector y at line %d\n", lineNum);
331     exit(8);
332     }
333    
334     dumpArray[j].v[i].y = atof( foo );
335    
336     foo = strtok(NULL, " ,;\t");
337     if(foo == NULL){
338     printf("error in reading vector z at line %d\n", lineNum);
339     exit(8);
340     }
341    
342     dumpArray[j].v[i].z = atof( foo );
343     }
344     }
345    
346 mmeineke 38 }
347     }
348    
349     (void)fclose(in_file);
350    
351     printf( "done\n"
352     "\n" );
353     fflush( stdout );
354    
355    
356    
357     // do calculations here.
358    
359 mmeineke 43 if( calcGofR ){
360    
361     fprintf( stdout,
362     "Calculating the g(r) between atoms \"%s\" and \"%s\"...",
363     gofR1, gofR2 );
364     fflush( stdout );
365    
366     // gofr call
367 mmeineke 82 GofR( out_prefix, gofR1, gofR2, dumpArray, nFrames, startFrame, endFrame );
368 mmeineke 43
369     fprintf( stdout,
370     " done.\n"
371     "\n");
372     fflush(stdout);
373     }
374 mmeineke 38
375 mmeineke 43 if( calcRMSD ){
376    
377     fprintf( stdout,
378     "Calculating the RMSD..." );
379     fflush( stdout );
380    
381     // RMSD call
382 mmeineke 38
383 mmeineke 43
384     fprintf( stdout,
385     " done.\n"
386     "\n");
387     fflush(stdout);
388     }
389 mmeineke 38
390 mmeineke 43 if( calcMuCorr ){
391    
392     fprintf( stdout,
393     "Calculating the mu correlation for \"%s\"...",
394     muCorr);
395     fflush( stdout );
396    
397     // muCorr call
398 mmeineke 38
399 mmeineke 43
400     fprintf( stdout,
401     " done.\n"
402     "\n");
403     fflush(stdout);
404     }
405    
406     if( calcCosCorr ){
407    
408     fprintf( stdout,
409     "Calculating the cos correlation between \"%s\" and \"%s\"...",
410     cosCorr1, cosCorr2 );
411     fflush( stdout );
412    
413 mmeineke 82 cosCorr( out_prefix, cosCorr1, cosCorr2, dumpArray, nFrames, startFrame,
414     endFrame );
415 mmeineke 43
416     fprintf( stdout,
417     " done.\n"
418     "\n");
419     fflush(stdout);
420     }
421 mmeineke 38
422     return 0;
423    
424     }
425    
426    
427 mmeineke 46 void map( double *x, double *y, double *z,
428     double boxX, double boxY, double boxZ ){
429    
430     *x -= boxX * copysign(1.0,*x) * floor( fabs( *x/boxX ) + 0.5 );
431     *y -= boxY * copysign(1.0,*y) * floor( fabs( *y/boxY ) + 0.5 );
432     *z -= boxZ * copysign(1.0,*z) * floor( fabs( *z/boxZ ) + 0.5 );
433 mmeineke 38
434 mmeineke 46 }
435    
436    
437 mmeineke 38 /***************************************************************************
438     * prints out the usage for the command line arguments, then exits.
439     ***************************************************************************/
440    
441     void usage(){
442 mmeineke 43 (void)fprintf(stdout,
443     "The proper usage is: %s [options] <xyz_file>\n"
444     "\n"
445     "Options:\n"
446     "\n"
447     " short:\n"
448     " ------\n"
449     " -h Display this message\n"
450     " -o <prefix> The output prefix\n"
451     " -r Calculate the RMSD\n"
452     " -g Calculate all to all g(r)\n"
453    
454     "\n"
455     " long:\n"
456     " -----\n"
457     " --GofR <atom1> <atom2> Calculates g(r) between atom1 and atom 2\n"
458     " -note: \"all\" will do all atoms\n"
459     " --MuCorr <atom> Calculate mu correlation of atom\n"
460     " --CosCorr <atom1> <atom2> Calculate the cos correlation between atom1 and atom2\n"
461 mmeineke 82 " --startFrame <frame#> Specifies a frame to start correlating\n"
462     " --endFrame <frame#> Specifies a frame to stop correlating.\n"
463 mmeineke 43
464     "\n"
465     "\n",
466 mmeineke 38 program_name);
467     exit(8);
468     }