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, 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"
53 static int nFloatSize(gmx_trr_header_t* sh)
59 nflsize = sh->box_size / (DIM * DIM);
63 nflsize = sh->x_size / (sh->natoms * DIM);
67 nflsize = sh->v_size / (sh->natoms * DIM);
71 nflsize = sh->f_size / (sh->natoms * DIM);
75 gmx_file("Can not determine precision of trr file");
78 if (((nflsize != sizeof(float)) && (nflsize != sizeof(double))))
80 gmx_fatal(FARGS, "Float size %d. Maybe different CPU?", nflsize);
86 /* Returns whether a valid frame header was read. Upon exit, *bOK is
87 TRUE if a normal outcome resulted. Usually that is the same thing,
88 but encountering the end of the file before reading the magic
89 integer is a normal outcome for TRR reading, and it does not
90 produce a valid frame header, so the values differ in that case.
91 That does not exclude the possibility of a reading error between
92 frames, but the trajectory-handling infrastructure needs an
93 overhaul before we can handle that. */
94 static gmx_bool do_trr_frame_header(t_fileio* fio, bool bRead, gmx_trr_header_t* sh, gmx_bool* bOK)
96 const int magicValue = 1993;
97 int magic = magicValue;
98 static gmx_bool bFirst = TRUE;
103 if (!gmx_fio_do_int(fio, magic))
105 /* Failed to read an integer, which should be the magic
106 number, which usually means we've reached the end
107 of the file (but could be an I/O error that we now
111 if (magic != magicValue)
115 "Failed to find GROMACS magic number in trr frame header, so this is not a trr "
121 *bOK = *bOK && gmx_fio_do_string(fio, buf);
124 fprintf(stderr, "trr version: %s ", buf);
129 sprintf(buf, "GMX_trn_file");
130 *bOK = *bOK && gmx_fio_do_string(fio, buf);
132 *bOK = *bOK && gmx_fio_do_int(fio, sh->ir_size);
133 *bOK = *bOK && gmx_fio_do_int(fio, sh->e_size);
134 *bOK = *bOK && gmx_fio_do_int(fio, sh->box_size);
135 *bOK = *bOK && gmx_fio_do_int(fio, sh->vir_size);
136 *bOK = *bOK && gmx_fio_do_int(fio, sh->pres_size);
137 *bOK = *bOK && gmx_fio_do_int(fio, sh->top_size);
138 *bOK = *bOK && gmx_fio_do_int(fio, sh->sym_size);
139 *bOK = *bOK && gmx_fio_do_int(fio, sh->x_size);
140 *bOK = *bOK && gmx_fio_do_int(fio, sh->v_size);
141 *bOK = *bOK && gmx_fio_do_int(fio, sh->f_size);
142 *bOK = *bOK && gmx_fio_do_int(fio, sh->natoms);
148 sh->bDouble = (nFloatSize(sh) == sizeof(double));
149 gmx_fio_setprecision(fio, sh->bDouble);
153 fprintf(stderr, "(%s precision)\n", sh->bDouble ? "double" : "single");
157 /* Note that TRR wasn't defined to be extensible, so we can't fix
158 * the fact that we used a default int for the step number, which
159 * is typically defined to be signed and 32 bit. */
160 int intStep = sh->step;
161 *bOK = *bOK && gmx_fio_do_int(fio, intStep);
163 *bOK = *bOK && gmx_fio_do_int(fio, sh->nre);
164 *bOK = *bOK && gmx_fio_do_real(fio, sh->t);
165 *bOK = *bOK && gmx_fio_do_real(fio, sh->lambda);
170 static gmx_bool do_trr_frame_data(t_fileio* fio, gmx_trr_header_t* sh, rvec* box, rvec* x, rvec* v, rvec* f)
176 if (sh->box_size != 0)
178 bOK = bOK && gmx_fio_ndo_rvec(fio, box, DIM);
180 if (sh->vir_size != 0)
182 bOK = bOK && gmx_fio_ndo_rvec(fio, pv, DIM);
184 if (sh->pres_size != 0)
186 bOK = bOK && gmx_fio_ndo_rvec(fio, pv, DIM);
190 bOK = bOK && gmx_fio_ndo_rvec(fio, x, sh->natoms);
194 bOK = bOK && gmx_fio_ndo_rvec(fio, v, sh->natoms);
198 bOK = bOK && gmx_fio_ndo_rvec(fio, f, sh->natoms);
204 static gmx_bool do_trr_frame(t_fileio* fio,
215 gmx_trr_header_t* sh;
221 sh->box_size = (box) ? sizeof(matrix) : 0;
222 sh->x_size = ((x) ? (*natoms * sizeof(x[0])) : 0);
223 sh->v_size = ((v) ? (*natoms * sizeof(v[0])) : 0);
224 sh->f_size = ((f) ? (*natoms * sizeof(f[0])) : 0);
225 sh->natoms = *natoms;
229 sh->lambda = *lambda;
231 if (!do_trr_frame_header(fio, bRead, sh, &bOK))
237 *natoms = sh->natoms;
240 *lambda = sh->lambda;
243 gmx_file("inputrec in trr file");
247 gmx_file("energies in trr file");
251 gmx_file("topology in trr file");
255 gmx_file("symbol table in trr file");
258 bOK = do_trr_frame_data(fio, sh, box, x, v, f);
265 /************************************************************
267 * The following routines are the exported ones
269 ************************************************************/
271 void gmx_trr_read_single_header(const char* fn, gmx_trr_header_t* header)
273 t_fileio* fio = gmx_trr_open(fn, "r");
275 if (!do_trr_frame_header(fio, true, header, &bOK))
277 gmx_fatal(FARGS, "Empty file %s", fn);
282 gmx_bool gmx_trr_read_frame_header(t_fileio* fio, gmx_trr_header_t* header, gmx_bool* bOK)
284 return do_trr_frame_header(fio, true, header, bOK);
287 void gmx_trr_write_single_frame(const char* fn,
297 t_fileio* fio = gmx_trr_open(fn, "w");
298 do_trr_frame(fio, false, &step, &t, &lambda, const_cast<rvec*>(box), &natoms,
299 const_cast<rvec*>(x), const_cast<rvec*>(v), const_cast<rvec*>(f));
303 void gmx_trr_read_single_frame(const char* fn,
313 t_fileio* fio = gmx_trr_open(fn, "r");
314 do_trr_frame(fio, true, step, t, lambda, box, natoms, x, v, f);
318 void gmx_trr_write_frame(t_fileio* fio,
328 if (!do_trr_frame(fio, false, &step, &t, &lambda, const_cast<rvec*>(box), &natoms,
329 const_cast<rvec*>(x), const_cast<rvec*>(v), const_cast<rvec*>(f)))
331 gmx_file("Cannot write trajectory frame; maybe you are out of disk space?");
336 gmx_bool gmx_trr_read_frame(t_fileio* fio,
346 return do_trr_frame(fio, true, step, t, lambda, box, natoms, x, v, f);
349 gmx_bool gmx_trr_read_frame_data(t_fileio* fio, gmx_trr_header_t* header, rvec* box, rvec* x, rvec* v, rvec* f)
351 return do_trr_frame_data(fio, header, box, x, v, f);
354 t_fileio* gmx_trr_open(const char* fn, const char* mode)
356 return gmx_fio_open(fn, mode);
359 void gmx_trr_close(t_fileio* fio)