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, 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"
51 #define GROMACS_MAGIC 1993
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
);
87 do_trr_frame_header(t_fileio
*fio
, bool bRead
, gmx_trr_header_t
*sh
, gmx_bool
*bOK
)
89 int magic
= GROMACS_MAGIC
;
90 static gmx_bool bFirst
= TRUE
;
95 if (!gmx_fio_do_int(fio
, magic
) || magic
!= GROMACS_MAGIC
)
102 *bOK
= *bOK
&& gmx_fio_do_string(fio
, buf
);
105 fprintf(stderr
, "trr 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
);
147 do_trr_frame_data(t_fileio
*fio
, gmx_trr_header_t
*sh
,
148 rvec
*box
, rvec
*x
, rvec
*v
, rvec
*f
)
154 if (sh
->box_size
!= 0)
156 bOK
= bOK
&& gmx_fio_ndo_rvec(fio
, box
, DIM
);
158 if (sh
->vir_size
!= 0)
160 bOK
= bOK
&& gmx_fio_ndo_rvec(fio
, pv
, DIM
);
162 if (sh
->pres_size
!= 0)
164 bOK
= bOK
&& gmx_fio_ndo_rvec(fio
, pv
, DIM
);
168 bOK
= bOK
&& gmx_fio_ndo_rvec(fio
, x
, sh
->natoms
);
172 bOK
= bOK
&& gmx_fio_ndo_rvec(fio
, v
, sh
->natoms
);
176 bOK
= bOK
&& gmx_fio_ndo_rvec(fio
, f
, sh
->natoms
);
183 do_trr_frame(t_fileio
*fio
, bool bRead
, int *step
, real
*t
, real
*lambda
,
184 rvec
*box
, int *natoms
, rvec
*x
, rvec
*v
, rvec
*f
)
186 gmx_trr_header_t
*sh
;
192 sh
->box_size
= (box
) ? sizeof(matrix
) : 0;
193 sh
->x_size
= ((x
) ? (*natoms
*sizeof(x
[0])) : 0);
194 sh
->v_size
= ((v
) ? (*natoms
*sizeof(v
[0])) : 0);
195 sh
->f_size
= ((f
) ? (*natoms
*sizeof(f
[0])) : 0);
196 sh
->natoms
= *natoms
;
200 sh
->lambda
= *lambda
;
202 if (!do_trr_frame_header(fio
, bRead
, sh
, &bOK
))
208 *natoms
= sh
->natoms
;
211 *lambda
= sh
->lambda
;
214 gmx_file("inputrec in trr file");
218 gmx_file("energies in trr file");
222 gmx_file("topology in trr file");
226 gmx_file("symbol table in trr file");
229 bOK
= do_trr_frame_data(fio
, sh
, box
, x
, v
, f
);
236 /************************************************************
238 * The following routines are the exported ones
240 ************************************************************/
242 void gmx_trr_read_single_header(const char *fn
, gmx_trr_header_t
*header
)
244 t_fileio
*fio
= gmx_trr_open(fn
, "r");
246 if (!do_trr_frame_header(fio
, true, header
, &bOK
))
248 gmx_fatal(FARGS
, "Empty file %s", fn
);
253 gmx_bool
gmx_trr_read_frame_header(t_fileio
*fio
, gmx_trr_header_t
*header
, gmx_bool
*bOK
)
255 return do_trr_frame_header(fio
, true, header
, bOK
);
258 void gmx_trr_write_single_frame(const char *fn
, int step
, real t
, real lambda
,
259 rvec
*box
, int natoms
, rvec
*x
, rvec
*v
, rvec
*f
)
261 t_fileio
*fio
= gmx_trr_open(fn
, "w");
262 do_trr_frame(fio
, false, &step
, &t
, &lambda
, box
, &natoms
, x
, v
, f
);
266 void gmx_trr_read_single_frame(const char *fn
, int *step
, real
*t
, real
*lambda
,
267 rvec
*box
, int *natoms
, rvec
*x
, rvec
*v
, rvec
*f
)
269 t_fileio
*fio
= gmx_trr_open(fn
, "r");
270 do_trr_frame(fio
, true, step
, t
, lambda
, box
, natoms
, x
, v
, f
);
274 void gmx_trr_write_frame(t_fileio
*fio
, int step
, real t
, real lambda
,
275 rvec
*box
, int natoms
, rvec
*x
, rvec
*v
, rvec
*f
)
277 if (!do_trr_frame(fio
, false, &step
, &t
, &lambda
, box
, &natoms
, x
, v
, f
))
279 gmx_file("Cannot write trajectory frame; maybe you are out of disk space?");
284 gmx_bool
gmx_trr_read_frame(t_fileio
*fio
, int *step
, real
*t
, real
*lambda
,
285 rvec
*box
, int *natoms
, rvec
*x
, rvec
*v
, rvec
*f
)
287 return do_trr_frame(fio
, true, step
, t
, lambda
, box
, natoms
, x
, v
, f
);
290 gmx_bool
gmx_trr_read_frame_data(t_fileio
*fio
, gmx_trr_header_t
*header
,
291 rvec
*box
, rvec
*x
, rvec
*v
, rvec
*f
)
293 return do_trr_frame_data(fio
, header
, box
, x
, v
, f
);
296 t_fileio
*gmx_trr_open(const char *fn
, const char *mode
)
298 return gmx_fio_open(fn
, mode
);
301 void gmx_trr_close(t_fileio
*fio
)