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, by the GROMACS development team, led by
7 * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
8 * and including many others, as listed in the AUTHORS file in the
9 * top-level source directory and at http://www.gromacs.org.
11 * GROMACS is free software; you can redistribute it and/or
12 * modify it under the terms of the GNU Lesser General Public License
13 * as published by the Free Software Foundation; either version 2.1
14 * of the License, or (at your option) any later version.
16 * GROMACS is distributed in the hope that it will be useful,
17 * but WITHOUT ANY WARRANTY; without even the implied warranty of
18 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
19 * Lesser General Public License for more details.
21 * You should have received a copy of the GNU Lesser General Public
22 * License along with GROMACS; if not, see
23 * http://www.gnu.org/licenses, or write to the Free Software Foundation,
24 * Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA.
26 * If you want to redistribute modifications to GROMACS, please
27 * consider that scientific software is very special. Version
28 * control is crucial - bugs must be traceable. We will be happy to
29 * consider code for inclusion in the official distribution, but
30 * derived work must not be called official GROMACS. Details are found
31 * in the README & COPYING files - if they are missing, get the
32 * official version at http://www.gromacs.org.
34 * To help us fund GROMACS development, we humbly ask that you cite
35 * the research papers on the package. Check out http://www.gromacs.org.
43 #include "gromacs/fileio/gmxfio.h"
44 #include "gromacs/fileio/gmxfio-xdr.h"
45 #include "gromacs/mdtypes/md_enums.h"
46 #include "gromacs/utility/fatalerror.h"
47 #include "gromacs/utility/futil.h"
48 #include "gromacs/utility/smalloc.h"
52 static int nFloatSize(gmx_trr_header_t
*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 trr file");
77 if (((nflsize
!= sizeof(float)) && (nflsize
!= sizeof(double))))
79 gmx_fatal(FARGS
, "Float size %d. Maybe different CPU?", nflsize
);
85 /* Returns whether a valid frame header was read. Upon exit, *bOK is
86 TRUE if a normal outcome resulted. Usually that is the same thing,
87 but encountering the end of the file before reading the magic
88 integer is a normal outcome for TRR reading, and it does not
89 produce a valid frame header, so the values differ in that case.
90 That does not exclude the possibility of a reading error between
91 frames, but the trajectory-handling infrastructure needs an
92 overhaul before we can handle that. */
94 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
)
114 gmx_fatal(FARGS
, "Failed to find GROMACS magic number in trr frame header, so this is not a trr file!\n");
120 *bOK
= *bOK
&& gmx_fio_do_string(fio
, buf
);
123 fprintf(stderr
, "trr version: %s ", buf
);
128 sprintf(buf
, "GMX_trn_file");
129 *bOK
= *bOK
&& gmx_fio_do_string(fio
, buf
);
131 *bOK
= *bOK
&& gmx_fio_do_int(fio
, sh
->ir_size
);
132 *bOK
= *bOK
&& gmx_fio_do_int(fio
, sh
->e_size
);
133 *bOK
= *bOK
&& gmx_fio_do_int(fio
, sh
->box_size
);
134 *bOK
= *bOK
&& gmx_fio_do_int(fio
, sh
->vir_size
);
135 *bOK
= *bOK
&& gmx_fio_do_int(fio
, sh
->pres_size
);
136 *bOK
= *bOK
&& gmx_fio_do_int(fio
, sh
->top_size
);
137 *bOK
= *bOK
&& gmx_fio_do_int(fio
, sh
->sym_size
);
138 *bOK
= *bOK
&& gmx_fio_do_int(fio
, sh
->x_size
);
139 *bOK
= *bOK
&& gmx_fio_do_int(fio
, sh
->v_size
);
140 *bOK
= *bOK
&& gmx_fio_do_int(fio
, sh
->f_size
);
141 *bOK
= *bOK
&& gmx_fio_do_int(fio
, sh
->natoms
);
147 sh
->bDouble
= (nFloatSize(sh
) == sizeof(double));
148 gmx_fio_setprecision(fio
, sh
->bDouble
);
152 fprintf(stderr
, "(%s precision)\n", sh
->bDouble
? "double" : "single");
156 /* Note that TRR wasn't defined to be extensible, so we can't fix
157 * the fact that we used a default int for the step number, which
158 * is typically defined to be signed and 32 bit. */
159 int intStep
= sh
->step
;
160 *bOK
= *bOK
&& gmx_fio_do_int(fio
, intStep
);
162 *bOK
= *bOK
&& gmx_fio_do_int(fio
, sh
->nre
);
163 *bOK
= *bOK
&& gmx_fio_do_real(fio
, sh
->t
);
164 *bOK
= *bOK
&& gmx_fio_do_real(fio
, sh
->lambda
);
170 do_trr_frame_data(t_fileio
*fio
, gmx_trr_header_t
*sh
,
171 rvec
*box
, rvec
*x
, rvec
*v
, rvec
*f
)
177 if (sh
->box_size
!= 0)
179 bOK
= bOK
&& gmx_fio_ndo_rvec(fio
, box
, DIM
);
181 if (sh
->vir_size
!= 0)
183 bOK
= bOK
&& gmx_fio_ndo_rvec(fio
, pv
, DIM
);
185 if (sh
->pres_size
!= 0)
187 bOK
= bOK
&& gmx_fio_ndo_rvec(fio
, pv
, DIM
);
191 bOK
= bOK
&& gmx_fio_ndo_rvec(fio
, x
, sh
->natoms
);
195 bOK
= bOK
&& gmx_fio_ndo_rvec(fio
, v
, sh
->natoms
);
199 bOK
= bOK
&& gmx_fio_ndo_rvec(fio
, f
, sh
->natoms
);
206 do_trr_frame(t_fileio
*fio
, bool bRead
, gmx_int64_t
*step
, real
*t
, real
*lambda
,
207 rvec
*box
, int *natoms
, rvec
*x
, rvec
*v
, rvec
*f
)
209 gmx_trr_header_t
*sh
;
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_trr_frame_header(fio
, bRead
, sh
, &bOK
))
231 *natoms
= sh
->natoms
;
234 *lambda
= sh
->lambda
;
237 gmx_file("inputrec in trr file");
241 gmx_file("energies in trr file");
245 gmx_file("topology in trr file");
249 gmx_file("symbol table in trr file");
252 bOK
= do_trr_frame_data(fio
, sh
, box
, x
, v
, f
);
259 /************************************************************
261 * The following routines are the exported ones
263 ************************************************************/
265 void gmx_trr_read_single_header(const char *fn
, gmx_trr_header_t
*header
)
267 t_fileio
*fio
= gmx_trr_open(fn
, "r");
269 if (!do_trr_frame_header(fio
, true, header
, &bOK
))
271 gmx_fatal(FARGS
, "Empty file %s", fn
);
276 gmx_bool
gmx_trr_read_frame_header(t_fileio
*fio
, gmx_trr_header_t
*header
, gmx_bool
*bOK
)
278 return do_trr_frame_header(fio
, true, header
, bOK
);
281 void gmx_trr_write_single_frame(const char *fn
, gmx_int64_t step
, real t
, real lambda
,
282 const rvec
*box
, int natoms
, const rvec
*x
, const rvec
*v
, const rvec
*f
)
284 t_fileio
*fio
= gmx_trr_open(fn
, "w");
285 do_trr_frame(fio
, false, &step
, &t
, &lambda
, const_cast<rvec
*>(box
), &natoms
, const_cast<rvec
*>(x
), const_cast<rvec
*>(v
), const_cast<rvec
*>(f
));
289 void gmx_trr_read_single_frame(const char *fn
, gmx_int64_t
*step
, real
*t
, real
*lambda
,
290 rvec
*box
, int *natoms
, rvec
*x
, rvec
*v
, rvec
*f
)
292 t_fileio
*fio
= gmx_trr_open(fn
, "r");
293 do_trr_frame(fio
, true, step
, t
, lambda
, box
, natoms
, x
, v
, f
);
297 void gmx_trr_write_frame(t_fileio
*fio
, gmx_int64_t step
, real t
, real lambda
,
298 const rvec
*box
, int natoms
, const rvec
*x
, const rvec
*v
, const rvec
*f
)
300 if (!do_trr_frame(fio
, false, &step
, &t
, &lambda
, const_cast<rvec
*>(box
), &natoms
, const_cast<rvec
*>(x
), const_cast<rvec
*>(v
), const_cast<rvec
*>(f
)))
302 gmx_file("Cannot write trajectory frame; maybe you are out of disk space?");
307 gmx_bool
gmx_trr_read_frame(t_fileio
*fio
, gmx_int64_t
*step
, real
*t
, real
*lambda
,
308 rvec
*box
, int *natoms
, rvec
*x
, rvec
*v
, rvec
*f
)
310 return do_trr_frame(fio
, true, step
, t
, lambda
, box
, natoms
, x
, v
, f
);
313 gmx_bool
gmx_trr_read_frame_data(t_fileio
*fio
, gmx_trr_header_t
*header
,
314 rvec
*box
, rvec
*x
, rvec
*v
, rvec
*f
)
316 return do_trr_frame_data(fio
, header
, box
, x
, v
, f
);
319 t_fileio
*gmx_trr_open(const char *fn
, const char *mode
)
321 return gmx_fio_open(fn
, mode
);
324 void gmx_trr_close(t_fileio
*fio
)