2 * This file is part of the GROMACS molecular simulation package.
4 * Copyright (c) 1991-2000, University of Groningen, The Netherlands.
5 * Copyright (c) 2001-2004, The GROMACS development team.
6 * Copyright (c) 2013,2014,2015,2016,2018 by the GROMACS development team.
7 * Copyright (c) 2019,2020,2021, by the GROMACS development team, led by
8 * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
9 * and including many others, as listed in the AUTHORS file in the
10 * top-level source directory and at http://www.gromacs.org.
12 * GROMACS is free software; you can redistribute it and/or
13 * modify it under the terms of the GNU Lesser General Public License
14 * as published by the Free Software Foundation; either version 2.1
15 * of the License, or (at your option) any later version.
17 * GROMACS is distributed in the hope that it will be useful,
18 * but WITHOUT ANY WARRANTY; without even the implied warranty of
19 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
20 * Lesser General Public License for more details.
22 * You should have received a copy of the GNU Lesser General Public
23 * License along with GROMACS; if not, see
24 * http://www.gnu.org/licenses, or write to the Free Software Foundation,
25 * Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA.
27 * If you want to redistribute modifications to GROMACS, please
28 * consider that scientific software is very special. Version
29 * control is crucial - bugs must be traceable. We will be happy to
30 * consider code for inclusion in the official distribution, but
31 * derived work must not be called official GROMACS. Details are found
32 * in the README & COPYING files - if they are missing, get the
33 * official version at http://www.gromacs.org.
35 * To help us fund GROMACS development, we humbly ask that you cite
36 * the research papers on the package. Check out http://www.gromacs.org.
44 #include "gromacs/fileio/gmxfio.h"
45 #include "gromacs/fileio/gmxfio_xdr.h"
46 #include "gromacs/mdtypes/md_enums.h"
47 #include "gromacs/utility/fatalerror.h"
48 #include "gromacs/utility/futil.h"
49 #include "gromacs/utility/smalloc.h"
51 static int nFloatSize(gmx_trr_header_t* sh)
57 nflsize = sh->box_size / (DIM * DIM);
61 nflsize = sh->x_size / (sh->natoms * DIM);
65 nflsize = sh->v_size / (sh->natoms * DIM);
69 nflsize = sh->f_size / (sh->natoms * DIM);
73 gmx_file("Can not determine precision of trr file");
76 if (((nflsize != sizeof(float)) && (nflsize != sizeof(double))))
78 gmx_fatal(FARGS, "Float size %d. Maybe different CPU?", nflsize);
84 /* Returns whether a valid frame header was read. Upon exit, *bOK is
85 TRUE if a normal outcome resulted. Usually that is the same thing,
86 but encountering the end of the file before reading the magic
87 integer is a normal outcome for TRR reading, and it does not
88 produce a valid frame header, so the values differ in that case.
89 That does not exclude the possibility of a reading error between
90 frames, but the trajectory-handling infrastructure needs an
91 overhaul before we can handle that. */
92 static gmx_bool do_trr_frame_header(t_fileio* fio, bool bRead, gmx_trr_header_t* sh, gmx_bool* bOK)
94 const int magicValue = 1993;
95 int magic = magicValue;
96 static gmx_bool bFirst = TRUE;
101 if (!gmx_fio_do_int(fio, magic))
103 /* Failed to read an integer, which should be the magic
104 number, which usually means we've reached the end
105 of the file (but could be an I/O error that we now
109 if (magic != magicValue)
113 "Failed to find GROMACS magic number in trr frame header, so this is not a trr "
119 *bOK = *bOK && gmx_fio_do_string(fio, buf);
122 fprintf(stderr, "trr version: %s ", buf);
127 sprintf(buf, "GMX_trn_file");
128 *bOK = *bOK && gmx_fio_do_string(fio, buf);
130 *bOK = *bOK && gmx_fio_do_int(fio, sh->ir_size);
131 *bOK = *bOK && gmx_fio_do_int(fio, sh->e_size);
132 *bOK = *bOK && gmx_fio_do_int(fio, sh->box_size);
133 *bOK = *bOK && gmx_fio_do_int(fio, sh->vir_size);
134 *bOK = *bOK && gmx_fio_do_int(fio, sh->pres_size);
135 *bOK = *bOK && gmx_fio_do_int(fio, sh->top_size);
136 *bOK = *bOK && gmx_fio_do_int(fio, sh->sym_size);
137 *bOK = *bOK && gmx_fio_do_int(fio, sh->x_size);
138 *bOK = *bOK && gmx_fio_do_int(fio, sh->v_size);
139 *bOK = *bOK && gmx_fio_do_int(fio, sh->f_size);
140 *bOK = *bOK && gmx_fio_do_int(fio, sh->natoms);
146 sh->bDouble = (nFloatSize(sh) == sizeof(double));
147 gmx_fio_setprecision(fio, sh->bDouble);
151 fprintf(stderr, "(%s precision)\n", sh->bDouble ? "double" : "single");
155 /* Note that TRR wasn't defined to be extensible, so we can't fix
156 * the fact that we used a default int for the step number, which
157 * is typically defined to be signed and 32 bit. */
158 int intStep = sh->step;
159 *bOK = *bOK && gmx_fio_do_int(fio, intStep);
161 *bOK = *bOK && gmx_fio_do_int(fio, sh->nre);
162 *bOK = *bOK && gmx_fio_do_real(fio, sh->t);
163 *bOK = *bOK && gmx_fio_do_real(fio, sh->lambda);
168 static gmx_bool do_trr_frame_data(t_fileio* fio, gmx_trr_header_t* sh, rvec* box, rvec* x, rvec* v, rvec* f)
174 if (sh->box_size != 0)
176 bOK = bOK && gmx_fio_ndo_rvec(fio, box, DIM);
178 if (sh->vir_size != 0)
180 bOK = bOK && gmx_fio_ndo_rvec(fio, pv, DIM);
182 if (sh->pres_size != 0)
184 bOK = bOK && gmx_fio_ndo_rvec(fio, pv, DIM);
188 bOK = bOK && gmx_fio_ndo_rvec(fio, x, sh->natoms);
192 bOK = bOK && gmx_fio_ndo_rvec(fio, v, sh->natoms);
196 bOK = bOK && gmx_fio_ndo_rvec(fio, f, sh->natoms);
202 static gmx_bool do_trr_frame(t_fileio* fio,
213 gmx_trr_header_t* sh;
219 sh->box_size = (box) ? sizeof(matrix) : 0;
220 sh->x_size = ((x) ? (*natoms * sizeof(x[0])) : 0);
221 sh->v_size = ((v) ? (*natoms * sizeof(v[0])) : 0);
222 sh->f_size = ((f) ? (*natoms * sizeof(f[0])) : 0);
223 sh->natoms = *natoms;
227 sh->lambda = *lambda;
229 if (!do_trr_frame_header(fio, bRead, sh, &bOK))
235 *natoms = sh->natoms;
238 *lambda = sh->lambda;
241 gmx_file("inputrec in trr file");
245 gmx_file("energies in trr file");
249 gmx_file("topology in trr file");
253 gmx_file("symbol table in trr file");
256 bOK = do_trr_frame_data(fio, sh, box, x, v, f);
263 /************************************************************
265 * The following routines are the exported ones
267 ************************************************************/
269 void gmx_trr_read_single_header(const char* fn, gmx_trr_header_t* header)
271 t_fileio* fio = gmx_trr_open(fn, "r");
273 if (!do_trr_frame_header(fio, true, header, &bOK))
275 gmx_fatal(FARGS, "Empty file %s", fn);
280 gmx_bool gmx_trr_read_frame_header(t_fileio* fio, gmx_trr_header_t* header, gmx_bool* bOK)
282 return do_trr_frame_header(fio, true, header, bOK);
285 void gmx_trr_write_single_frame(const char* fn,
295 t_fileio* fio = gmx_trr_open(fn, "w");
301 const_cast<rvec*>(box),
303 const_cast<rvec*>(x),
304 const_cast<rvec*>(v),
305 const_cast<rvec*>(f));
309 void gmx_trr_read_single_frame(const char* fn,
319 t_fileio* fio = gmx_trr_open(fn, "r");
320 do_trr_frame(fio, true, step, t, lambda, box, natoms, x, v, f);
324 void gmx_trr_write_frame(t_fileio* fio,
334 if (!do_trr_frame(fio,
339 const_cast<rvec*>(box),
341 const_cast<rvec*>(x),
342 const_cast<rvec*>(v),
343 const_cast<rvec*>(f)))
345 gmx_file("Cannot write trajectory frame; maybe you are out of disk space?");
350 gmx_bool gmx_trr_read_frame(t_fileio* fio,
360 return do_trr_frame(fio, true, step, t, lambda, box, natoms, x, v, f);
363 gmx_bool gmx_trr_read_frame_data(t_fileio* fio, gmx_trr_header_t* header, rvec* box, rvec* x, rvec* v, rvec* f)
365 return do_trr_frame_data(fio, header, box, x, v, f);
368 t_fileio* gmx_trr_open(const char* fn, const char* mode)
370 return gmx_fio_open(fn, mode);
373 void gmx_trr_close(t_fileio* fio)