Fix implicit declarations in C files
[namd.git] / src / sortreplicas.c
1
2 #include "largefiles.h"  /* must be first! */
3
4 #include "vmdplugin.h"
5
6 extern int molfile_dcdplugin_init(void);
7 extern int molfile_dcdplugin_register(void *, vmdplugin_register_cb);
8 extern int molfile_dcdplugin_fini(void);
9 /*
10 extern int molfile_jsplugin_init(void);
11 extern int molfile_jsplugin_register(void *, vmdplugin_register_cb);
12 extern int molfile_jsplugin_fini(void);
13 */
14
15 #include "molfile_plugin.h"
16
17 static int register_cb(void *v, vmdplugin_t *p) {
18   *((vmdplugin_t **)v) = p;
19   return VMDPLUGIN_SUCCESS;
20 }
21
22 #include <sys/types.h>
23 #include <stdio.h>
24 #include <stdlib.h>
25 #include <string.h>
26 #include <errno.h>
27
28 int main(int argc, char **argv) {
29
30   molfile_timestep_t frame;
31   molfile_plugin_t *plugin;
32   char *output_root;
33   char *filename;
34   int num_replicas;
35   int runs_per_frame;
36   long long int final_step = -1;
37   long long int checkstep = -1;
38   int colvars;
39   FILE **hist_in;
40   FILE **hist_out;
41   FILE **colv_in;
42   FILE **colv_out;
43   void **traj_in;
44   void **traj_out;
45   int natoms=MOLFILE_NUMATOMS_UNKNOWN;
46   int i, i_run;
47
48   molfile_dcdplugin_init();
49   molfile_dcdplugin_register(&plugin, register_cb);
50
51   if ( argc < 4 || argc > 5 ) {
52     fprintf(stderr, "args: <job_output_root> <num_replicas> <runs_per_frame> [final_step]\n");
53     exit(-1);
54   }
55   output_root = argv[1];
56   num_replicas = atoi(argv[2]);
57   runs_per_frame = atoi(argv[3]);
58   if ( argc > 4 ) {
59     sscanf(argv[4], "%lld", &final_step);
60   }
61
62   filename = (char*) malloc(strlen(output_root)+100);
63   hist_in = (FILE**) malloc(num_replicas*sizeof(FILE*));
64   hist_out = (FILE**) malloc(num_replicas*sizeof(FILE*));
65   colv_in = (FILE**) malloc(num_replicas*sizeof(FILE*));
66   colv_out = (FILE**) malloc(num_replicas*sizeof(FILE*));
67   traj_in = (void**) malloc(num_replicas*sizeof(FILE*));
68   traj_out = (void**) malloc(num_replicas*sizeof(FILE*));
69
70   for ( i=0; i<num_replicas; ++i ) {
71     char *root_end;
72     if ( strstr(output_root,"%s") ) {
73       char istr[10];
74       sprintf(istr,"%d",i);
75       sprintf(filename,output_root,istr);
76     } else {
77       sprintf(filename,output_root,i);
78     }
79     root_end = filename + strlen(filename);
80
81     sprintf(root_end,".%d.history",i);
82     hist_in[i] = fopen(filename,"r");
83     if ( ! hist_in[i] ) {
84       fprintf(stderr, "error opening input file %s: %s\n",
85                                         filename, strerror(errno));
86       exit(-1);
87     }
88     sprintf(root_end,".%d.dcd",i);
89     traj_in[i] = plugin->open_file_read(filename,"dcd",&natoms);
90     if ( ! traj_in[i] ) {
91       fprintf(stderr, "error opening input file %s: %s\n",
92                                         filename, strerror(errno));
93       exit(-1);
94     }
95     sprintf(root_end,".%d.colvars.traj",i);
96     colv_in[i] = fopen(filename,"r");
97     if ( colv_in[i] ) {
98       if ( i == 0 ) {
99         printf("Found first input colvars trajectory file %s.\n", filename);
100         colvars = 1;
101       } else if ( ! colvars ) {
102         fprintf(stderr, "missing input colvars trajectory files before %s\n", filename);
103         exit(-1);
104       }
105     } else {
106       if ( i == 0 ) {
107         colvars = 0;
108       } else if ( colvars ) {
109         fprintf(stderr, "error opening input colvars trajectory file %s: %s\n",
110                                         filename, strerror(errno));
111         exit(-1);
112       }
113     }
114   }
115
116   for ( i=0; i<num_replicas; ++i ) {
117     char *root_end;
118     if ( strstr(output_root,"%s") ) {
119       char istr[10];
120       sprintf(istr,"%d",i);
121       sprintf(filename,output_root,istr);
122     } else {
123       sprintf(filename,output_root,i);
124     }
125     root_end = filename + strlen(filename);
126
127     sprintf(root_end,".%d.sort.history",i);
128     hist_out[i] = fopen(filename,"w");
129     if ( ! hist_out[i] ) {
130       fprintf(stderr, "error opening output file %s: %s\n",
131                                         filename, strerror(errno));
132       exit(-1);
133     }
134     sprintf(root_end,".%d.sort.dcd",i);
135     traj_out[i] = plugin->open_file_write(filename,"dcd",natoms);
136     if ( ! traj_out[i] ) {
137       fprintf(stderr, "error opening output file %s: %s\n",
138                                         filename, strerror(errno));
139       exit(-1);
140     }
141     if ( colvars ) {
142       sprintf(root_end,".%d.sort.colvars.traj",i);
143       colv_out[i] = fopen(filename,"w");
144       if ( ! colv_out[i] ) {
145         fprintf(stderr, "error opening output file %s: %s\n",
146                                         filename, strerror(errno));
147         exit(-1);
148       }
149     }
150   }
151
152   frame.coords = (float*) malloc(3*natoms*sizeof(float));
153   frame.velocities = (float*) NULL;
154
155 #define LINE_MAX 10000
156
157   i_run = 0;
158   for ( ; 1; ++i_run ) { /* loop until read fails */
159     char line[LINE_MAX];
160     for ( i=0; i<num_replicas; ++i ) {
161       char *r;
162       char sav;
163       int f1,f2;
164       int rc;
165       int rep_id = -1;
166       long long int step;
167       r = fgets(line, LINE_MAX, hist_in[i]);
168       if ( ! r ) { break; }
169       rc = sscanf(line, "%lld %n%d%n", &step, &f1, &rep_id, &f2);
170       if ( rc != 2 ) {
171         fprintf(stderr,"Format error for replica %d at line %d: %s",
172                                                         i, i_run, line);
173         exit(-1);
174       }
175       if ( i == 0 ) {
176         if ( step <= checkstep ) {
177           fprintf(stderr,"Step out of order for replica %d at line %d: %s",
178                                                         i, i_run, line);
179           exit(-1);
180         }
181         checkstep = step;
182         if ( final_step >= 0 && checkstep > final_step ) {
183           printf("Stopping after final step %lld.\n", final_step);
184           break;
185         }
186       } else if ( step != checkstep ) {
187         fprintf(stderr,"Step mismatch for replica %d at line %d: %s",
188                                                         i, i_run, line);
189         exit(-1);
190       }
191       if ( rep_id < 0 || rep_id >= num_replicas ) {
192         fprintf(stderr,"Invalid replica ID for replica %d at line %d: %s",
193                                                         i, i_run, line);
194         exit(-1);
195       }
196       sav = line[f1];
197       line[f1] = 0;
198       fprintf(hist_out[rep_id],"%s%d%s",line,i,line+f2);
199       line[f1] = sav;
200       if ( colvars ) {
201        long long int oldcstep = -1;
202        while ( 1 ) {
203         long long int cstep;
204         char cline[LINE_MAX];
205 #ifdef WIN32
206         __int64 oldpos = _ftelli64(colv_in[i]);
207 #else
208         off_t oldpos = ftello(colv_in[i]);
209 #endif
210         r = fgets(cline, LINE_MAX, colv_in[i]);
211         if ( ! r ) { break; }
212         if ( cline[0] == '#' ) {
213           fprintf(colv_out[rep_id],"%s",cline);
214           continue;
215         } 
216         rc = sscanf(cline, "%lld", &cstep);
217         if ( rc != 1 ) {
218           fprintf(stderr,"Format error in colvar trajectory for replica %d: %s",
219                                                         i, cline);
220           exit(-1);
221         }
222         if ( cstep == oldcstep ) continue;  /* filter out repeats */
223         if ( cstep < oldcstep ) {
224           fprintf(stderr,"Step out of order in colvar trajectory for replica %d: %s",
225                                                         i, cline);
226           exit(-1);
227         }
228         if ( cstep > step ) {
229 #ifdef WIN32
230           _fseeki64(colv_in[i], oldpos, SEEK_SET);
231 #else
232           fseeko(colv_in[i], oldpos, SEEK_SET);
233 #endif
234           break;
235         }
236         if ( i_run != 0 || oldcstep != -1 ) {  /* skip first entry */
237           fprintf(colv_out[rep_id],"%s",cline);
238         }
239         oldcstep = cstep;
240        }
241       }
242       if ( (i_run+1) % runs_per_frame ) continue;
243       rc = plugin->read_next_timestep(traj_in[i],natoms,&frame);
244       if ( rc == MOLFILE_SUCCESS ) {
245         plugin->write_timestep(traj_out[rep_id],&frame);
246       } else {
247         fprintf(stderr,"Unable to read frame for replica %d at line %d: %s",
248                                                         i, i_run, line);
249         break;
250       }
251     }
252     if ( i < num_replicas ) {
253       printf("Processed %d runs.\n",i_run);
254       if ( i ) fprintf(stderr,"Uneven input lengths for replica %d at line %d: %s",
255                                                         i, i_run, line);
256       break;
257     }
258   }
259
260   free(frame.coords);
261
262   for ( i=0; i<num_replicas; ++i ) {
263     if ( fclose(hist_in[i]) ) {
264       fprintf(stderr, "error closing history input file %d: %s\n", i, strerror(errno));
265     }
266     plugin->close_file_read(traj_in[i]);
267     if ( fclose(hist_out[i]) ) {
268       fprintf(stderr, "error closing history output file %d: %s\n", i, strerror(errno));
269     }
270     plugin->close_file_write(traj_out[i]);
271     if ( colvars ) {
272       if ( fclose(colv_in[i]) ) {
273         fprintf(stderr, "error closing colvars input file %d: %s\n", i, strerror(errno));
274       }
275       if ( fclose(colv_out[i]) ) {
276         fprintf(stderr, "error closing colvars output file %d: %s\n", i, strerror(errno));
277       }
278     }
279   }
280
281   molfile_dcdplugin_fini();
282   exit(0);
283 }
284