3 * This source code is part of
7 * GROningen MAchine for Chemical Simulations
10 * Written by David van der Spoel, Erik Lindahl, Berk Hess, and others.
11 * Copyright (c) 1991-2000, University of Groningen, The Netherlands.
12 * Copyright (c) 2001-2004, The GROMACS development team,
13 * check out http://www.gromacs.org for more information.
15 * This program is free software; you can redistribute it and/or
16 * modify it under the terms of the GNU General Public License
17 * as published by the Free Software Foundation; either version 2
18 * of the License, or (at your option) any later version.
20 * If you want to redistribute modifications, please consider that
21 * scientific software is very special. Version control is crucial -
22 * bugs must be traceable. We will be happy to consider code for
23 * inclusion in the official distribution, but derived work must not
24 * be called official GROMACS. Details are found in the README & COPYING
25 * files - if they are missing, get the official version at www.gromacs.org.
27 * To help us fund GROMACS development, we humbly ask that you cite
28 * the papers on the package - you can find them in the top README file.
30 * For more info, check our website at http://www.gromacs.org
33 * GROningen Mixture of Alchemy and Childrens' Stories
42 #include "gmx_fatal.h"
50 #define GROMACS_MAGIC 1993
52 static int nFloatSize(t_trnheader *sh)
58 nflsize = sh->box_size/(DIM*DIM);
62 nflsize = sh->x_size/(sh->natoms*DIM);
66 nflsize = sh->v_size/(sh->natoms*DIM);
70 nflsize = sh->f_size/(sh->natoms*DIM);
74 gmx_file("Can not determine precision of trn file");
77 if (((nflsize != sizeof(float)) && (nflsize != sizeof(double))))
79 gmx_fatal(FARGS, "Float size %d. Maybe different CPU?", nflsize);
85 static gmx_bool do_trnheader(t_fileio *fio, gmx_bool bRead, t_trnheader *sh, gmx_bool *bOK)
87 int magic = GROMACS_MAGIC;
88 static gmx_bool bFirst = TRUE;
93 gmx_fio_checktype(fio);
95 if (!gmx_fio_do_int(fio, magic) || magic != GROMACS_MAGIC)
102 *bOK = *bOK && gmx_fio_do_string(fio, buf);
105 fprintf(stderr, "trn version: %s ", buf);
110 sprintf(buf, "GMX_trn_file");
111 *bOK = *bOK && gmx_fio_do_string(fio, buf);
113 *bOK = *bOK && gmx_fio_do_int(fio, sh->ir_size);
114 *bOK = *bOK && gmx_fio_do_int(fio, sh->e_size);
115 *bOK = *bOK && gmx_fio_do_int(fio, sh->box_size);
116 *bOK = *bOK && gmx_fio_do_int(fio, sh->vir_size);
117 *bOK = *bOK && gmx_fio_do_int(fio, sh->pres_size);
118 *bOK = *bOK && gmx_fio_do_int(fio, sh->top_size);
119 *bOK = *bOK && gmx_fio_do_int(fio, sh->sym_size);
120 *bOK = *bOK && gmx_fio_do_int(fio, sh->x_size);
121 *bOK = *bOK && gmx_fio_do_int(fio, sh->v_size);
122 *bOK = *bOK && gmx_fio_do_int(fio, sh->f_size);
123 *bOK = *bOK && gmx_fio_do_int(fio, sh->natoms);
129 sh->bDouble = (nFloatSize(sh) == sizeof(double));
130 gmx_fio_setprecision(fio, sh->bDouble);
134 fprintf(stderr, "(%s precision)\n", sh->bDouble ? "double" : "single");
138 *bOK = *bOK && gmx_fio_do_int(fio, sh->step);
139 *bOK = *bOK && gmx_fio_do_int(fio, sh->nre);
140 *bOK = *bOK && gmx_fio_do_real(fio, sh->t);
141 *bOK = *bOK && gmx_fio_do_real(fio, sh->lambda);
146 void pr_trnheader(FILE *fp, int indent, char *title, t_trnheader *sh)
150 indent = pr_title(fp, indent, title);
151 (void) pr_indent(fp, indent);
152 (void) fprintf(fp, "box_size = %d\n", sh->box_size);
153 (void) pr_indent(fp, indent);
154 (void) fprintf(fp, "x_size = %d\n", sh->x_size);
155 (void) pr_indent(fp, indent);
156 (void) fprintf(fp, "v_size = %d\n", sh->v_size);
157 (void) pr_indent(fp, indent);
158 (void) fprintf(fp, "f_size = %d\n", sh->f_size);
160 (void) pr_indent(fp, indent);
161 (void) fprintf(fp, "natoms = %d\n", sh->natoms);
162 (void) pr_indent(fp, indent);
163 (void) fprintf(fp, "step = %d\n", sh->step);
164 (void) pr_indent(fp, indent);
165 (void) fprintf(fp, "t = %e\n", sh->t);
166 (void) pr_indent(fp, indent);
167 (void) fprintf(fp, "lambda = %e\n", sh->lambda);
171 static gmx_bool do_htrn(t_fileio *fio, t_trnheader *sh,
172 rvec *box, rvec *x, rvec *v, rvec *f)
178 if (sh->box_size != 0)
180 bOK = bOK && gmx_fio_ndo_rvec(fio, box, DIM);
182 if (sh->vir_size != 0)
184 bOK = bOK && gmx_fio_ndo_rvec(fio, pv, DIM);
186 if (sh->pres_size != 0)
188 bOK = bOK && gmx_fio_ndo_rvec(fio, pv, DIM);
192 bOK = bOK && gmx_fio_ndo_rvec(fio, x, sh->natoms);
196 bOK = bOK && gmx_fio_ndo_rvec(fio, v, sh->natoms);
200 bOK = bOK && gmx_fio_ndo_rvec(fio, f, sh->natoms);
206 static gmx_bool do_trn(t_fileio *fio, gmx_bool bRead, int *step, real *t, real *lambda,
207 rvec *box, int *natoms, rvec *x, rvec *v, rvec *f)
215 sh->box_size = (box) ? sizeof(matrix) : 0;
216 sh->x_size = ((x) ? (*natoms*sizeof(x[0])) : 0);
217 sh->v_size = ((v) ? (*natoms*sizeof(v[0])) : 0);
218 sh->f_size = ((f) ? (*natoms*sizeof(f[0])) : 0);
219 sh->natoms = *natoms;
223 sh->lambda = *lambda;
225 if (!do_trnheader(fio, bRead, sh, &bOK))
231 *natoms = sh->natoms;
234 *lambda = sh->lambda;
237 gmx_file("inputrec in trn file");
241 gmx_file("energies in trn file");
245 gmx_file("topology in trn file");
249 gmx_file("symbol table in trn file");
252 bOK = do_htrn(fio, sh, box, x, v, f);
259 /************************************************************
261 * The following routines are the exported ones
263 ************************************************************/
265 void read_trnheader(const char *fn, t_trnheader *trn)
270 fio = open_trn(fn, "r");
271 if (!do_trnheader(fio, TRUE, trn, &bOK))
273 gmx_fatal(FARGS, "Empty file %s", fn);
278 gmx_bool fread_trnheader(t_fileio *fio, t_trnheader *trn, gmx_bool *bOK)
280 return do_trnheader(fio, TRUE, trn, bOK);
283 void write_trn(const char *fn, int step, real t, real lambda,
284 rvec *box, int natoms, rvec *x, rvec *v, rvec *f)
288 fio = open_trn(fn, "w");
289 do_trn(fio, FALSE, &step, &t, &lambda, box, &natoms, x, v, f);
293 void read_trn(const char *fn, int *step, real *t, real *lambda,
294 rvec *box, int *natoms, rvec *x, rvec *v, rvec *f)
298 fio = open_trn(fn, "r");
299 (void) do_trn(fio, TRUE, step, t, lambda, box, natoms, x, v, f);
303 void fwrite_trn(t_fileio *fio, int step, real t, real lambda,
304 rvec *box, int natoms, rvec *x, rvec *v, rvec *f)
306 if (do_trn(fio, FALSE, &step, &t, &lambda, box, &natoms, x, v, f) == FALSE)
308 gmx_file("Cannot write trajectory frame; maybe you are out of disk space?");
313 gmx_bool fread_trn(t_fileio *fio, int *step, real *t, real *lambda,
314 rvec *box, int *natoms, rvec *x, rvec *v, rvec *f)
316 return do_trn(fio, TRUE, step, t, lambda, box, natoms, x, v, f);
319 gmx_bool fread_htrn(t_fileio *fio, t_trnheader *trn, rvec *box, rvec *x, rvec *v,
322 return do_htrn(fio, trn, box, x, v, f);
325 t_fileio *open_trn(const char *fn, const char *mode)
327 return gmx_fio_open(fn, mode);
330 void close_trn(t_fileio *fio)