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,2017,2018,2019, 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.
47 #include "gromacs/fileio/checkpoint.h"
48 #include "gromacs/fileio/confio.h"
49 #include "gromacs/fileio/filetypes.h"
50 #include "gromacs/fileio/g96io.h"
51 #include "gromacs/fileio/gmxfio.h"
52 #include "gromacs/fileio/gmxfio_xdr.h"
53 #include "gromacs/fileio/groio.h"
54 #include "gromacs/fileio/oenv.h"
55 #include "gromacs/fileio/pdbio.h"
56 #include "gromacs/fileio/timecontrol.h"
57 #include "gromacs/fileio/tngio.h"
58 #include "gromacs/fileio/tpxio.h"
59 #include "gromacs/fileio/trrio.h"
60 #include "gromacs/fileio/xdrf.h"
61 #include "gromacs/fileio/xtcio.h"
62 #include "gromacs/math/vec.h"
63 #include "gromacs/mdtypes/md_enums.h"
64 #include "gromacs/topology/atoms.h"
65 #include "gromacs/topology/symtab.h"
66 #include "gromacs/topology/topology.h"
67 #include "gromacs/trajectory/trajectoryframe.h"
68 #include "gromacs/utility/fatalerror.h"
69 #include "gromacs/utility/futil.h"
70 #include "gromacs/utility/gmxassert.h"
71 #include "gromacs/utility/smalloc.h"
74 #include "gromacs/fileio/vmdio.h"
77 /* defines for frame counter output */
84 int flags
; /* flags for read_first/next_frame */
86 real t0
; /* time of the first frame, needed *
87 * for skipping frames with -dt */
88 real tf
; /* internal frame time */
91 gmx_tng_trajectory_t tng
;
95 char *persistent_line
; /* Persistent line for reading g96 trajectories */
97 gmx_vmdplugin_t
*vmdplugin
;
101 /* utility functions */
103 gmx_bool
bRmod_fd(double a
, double b
, double c
, gmx_bool bDouble
)
108 tol
= 2*(bDouble
? GMX_DOUBLE_EPS
: GMX_FLOAT_EPS
);
110 iq
= static_cast<int>((a
- b
+ tol
*a
)/c
);
112 return fabs(a
- b
- c
*iq
) <= tol
*fabs(a
);
117 int check_times2(real t
, real t0
, gmx_bool bDouble
)
122 /* since t is float, we can not use double precision for bRmod */
127 if ((!bTimeSet(TBEGIN
) || (t
>= rTimeValue(TBEGIN
))) &&
128 (!bTimeSet(TEND
) || (t
<= rTimeValue(TEND
))))
130 if (bTimeSet(TDELTA
) && !bRmod_fd(t
, t0
, rTimeValue(TDELTA
), bDouble
))
139 else if (bTimeSet(TEND
) && (t
>= rTimeValue(TEND
)))
145 fprintf(debug
, "t=%g, t0=%g, b=%g, e=%g, dt=%g: r=%d\n",
146 t
, t0
, rTimeValue(TBEGIN
), rTimeValue(TEND
), rTimeValue(TDELTA
), r
);
151 int check_times(real t
)
153 return check_times2(t
, t
, FALSE
);
156 static void initcount(t_trxstatus
*status
)
158 status
->__frame
= -1;
161 static void status_init(t_trxstatus
*status
)
164 status
->xframe
= nullptr;
165 status
->fio
= nullptr;
166 status
->__frame
= -1;
169 status
->persistent_line
= nullptr;
170 status
->tng
= nullptr;
174 int nframes_read(t_trxstatus
*status
)
176 return status
->__frame
;
179 static void printcount_(t_trxstatus
*status
, const gmx_output_env_t
*oenv
,
180 const char *l
, real t
)
182 if ((status
->__frame
< 2*SKIP1
|| status
->__frame
% SKIP1
== 0) &&
183 (status
->__frame
< 2*SKIP2
|| status
->__frame
% SKIP2
== 0) &&
184 (status
->__frame
< 2*SKIP3
|| status
->__frame
% SKIP3
== 0) &&
185 output_env_get_trajectory_io_verbosity(oenv
) != 0)
187 fprintf(stderr
, "\r%-14s %6d time %8.3f ", l
, status
->__frame
,
188 output_env_conv_time(oenv
, t
));
193 static void printcount(t_trxstatus
*status
, const gmx_output_env_t
*oenv
, real t
,
197 printcount_(status
, oenv
, bSkip
? "Skipping frame" : "Reading frame", t
);
200 static void printlast(t_trxstatus
*status
, const gmx_output_env_t
*oenv
, real t
)
202 printcount_(status
, oenv
, "Last frame", t
);
203 fprintf(stderr
, "\n");
207 static void printincomp(t_trxstatus
*status
, t_trxframe
*fr
)
209 if (fr
->not_ok
& HEADER_NOT_OK
)
211 fprintf(stderr
, "WARNING: Incomplete header: nr %d time %g\n",
212 status
->__frame
+1, fr
->time
);
216 fprintf(stderr
, "WARNING: Incomplete frame: nr %d time %g\n",
217 status
->__frame
+1, fr
->time
);
222 int prec2ndec(real prec
)
226 gmx_fatal(FARGS
, "DEATH HORROR prec (%g) <= 0 in prec2ndec", prec
);
229 return gmx::roundToInt(log(prec
)/log(10.0));
232 real
ndec2prec(int ndec
)
234 return pow(10.0, ndec
);
237 t_fileio
*trx_get_fileio(t_trxstatus
*status
)
242 float trx_get_time_of_final_frame(t_trxstatus
*status
)
244 t_fileio
*stfio
= trx_get_fileio(status
);
245 int filetype
= gmx_fio_getftp(stfio
);
249 if (filetype
== efXTC
)
252 xdr_xtc_get_last_frame_time(gmx_fio_getfp(stfio
),
253 gmx_fio_getxdr(stfio
),
254 status
->natoms
, &bOK
);
257 gmx_fatal(FARGS
, "Error reading last frame. Maybe seek not supported." );
260 else if (filetype
== efTNG
)
262 gmx_tng_trajectory_t tng
= status
->tng
;
265 gmx_fatal(FARGS
, "Error opening TNG file.");
267 lasttime
= gmx_tng_get_time_of_final_frame(tng
);
271 gmx_incons("Only supported for TNG and XTC");
276 void clear_trxframe(t_trxframe
*fr
, gmx_bool bFirst
)
282 fr
->bFepState
= FALSE
;
308 void set_trxframe_ePBC(t_trxframe
*fr
, int ePBC
)
310 fr
->bPBC
= (ePBC
== -1);
314 int write_trxframe_indexed(t_trxstatus
*status
, const t_trxframe
*fr
, int nind
,
315 const int *ind
, gmx_conect gc
)
318 rvec
*xout
= nullptr, *vout
= nullptr, *fout
= nullptr;
335 else if (status
->fio
)
337 ftp
= gmx_fio_getftp(status
->fio
);
341 gmx_incons("No input file available");
352 gmx_fatal(FARGS
, "Need coordinates to write a %s trajectory",
365 for (i
= 0; i
< nind
; i
++)
367 copy_rvec(fr
->v
[ind
[i
]], vout
[i
]);
373 for (i
= 0; i
< nind
; i
++)
375 copy_rvec(fr
->f
[ind
[i
]], fout
[i
]);
381 for (i
= 0; i
< nind
; i
++)
383 copy_rvec(fr
->x
[ind
[i
]], xout
[i
]);
391 for (i
= 0; i
< nind
; i
++)
393 copy_rvec(fr
->x
[ind
[i
]], xout
[i
]);
404 gmx_write_tng_from_trxframe(status
->tng
, fr
, nind
);
407 write_xtc(status
->fio
, nind
, fr
->step
, fr
->time
, fr
->box
, xout
, prec
);
410 gmx_trr_write_frame(status
->fio
, nframes_read(status
),
411 fr
->time
, fr
->step
, fr
->box
, nind
, xout
, vout
, fout
);
419 gmx_fatal(FARGS
, "Can not write a %s file without atom names",
422 sprintf(title
, "frame t= %.3f", fr
->time
);
425 write_hconf_indexed_p(gmx_fio_getfp(status
->fio
), title
, fr
->atoms
, nind
, ind
,
426 fr
->x
, fr
->bV
? fr
->v
: nullptr, fr
->box
);
430 write_pdbfile_indexed(gmx_fio_getfp(status
->fio
), title
, fr
->atoms
,
431 fr
->x
, -1, fr
->box
, ' ', fr
->step
, nind
, ind
, gc
, FALSE
);
435 sprintf(title
, "frame t= %.3f", fr
->time
);
436 write_g96_conf(gmx_fio_getfp(status
->fio
), title
, fr
, nind
, ind
);
439 gmx_fatal(FARGS
, "Sorry, write_trxframe_indexed can not write %s",
468 trjtools_gmx_prepare_tng_writing(const char *filename
,
473 const gmx_mtop_t
*mtop
,
474 gmx::ArrayRef
<const int> index
,
475 const char *index_group_name
)
477 if (filemode
!= 'w' && filemode
!= 'a')
479 gmx_incons("Sorry, can only prepare for TNG output.");
487 gmx_prepare_tng_writing(filename
,
496 else if ((infile
) && (efTNG
== fn2ftp(infile
)))
498 gmx_tng_trajectory_t tng_in
;
499 gmx_tng_open(infile
, 'r', &tng_in
);
501 gmx_prepare_tng_writing(filename
,
512 // we start from a file that is not a tng file or have been unable to load the
513 // input file, so we need to populate the fields independently of it
514 gmx_prepare_tng_writing(filename
,
526 void write_tng_frame(t_trxstatus
*status
,
529 gmx_write_tng_from_trxframe(status
->tng
, frame
, -1);
532 int write_trxframe(t_trxstatus
*status
, t_trxframe
*fr
, gmx_conect gc
)
548 gmx_tng_set_compression_precision(status
->tng
, prec
);
549 write_tng_frame(status
, fr
);
554 switch (gmx_fio_getftp(status
->fio
))
561 gmx_fatal(FARGS
, "Need coordinates to write a %s trajectory",
562 ftp2ext(gmx_fio_getftp(status
->fio
)));
567 switch (gmx_fio_getftp(status
->fio
))
570 write_xtc(status
->fio
, fr
->natoms
, fr
->step
, fr
->time
, fr
->box
, fr
->x
, prec
);
573 gmx_trr_write_frame(status
->fio
, fr
->step
, fr
->time
, fr
->lambda
, fr
->box
, fr
->natoms
,
574 fr
->bX
? fr
->x
: nullptr, fr
->bV
? fr
->v
: nullptr, fr
->bF
? fr
->f
: nullptr);
582 gmx_fatal(FARGS
, "Can not write a %s file without atom names",
583 ftp2ext(gmx_fio_getftp(status
->fio
)));
585 sprintf(title
, "frame t= %.3f", fr
->time
);
586 if (gmx_fio_getftp(status
->fio
) == efGRO
)
588 write_hconf_p(gmx_fio_getfp(status
->fio
), title
, fr
->atoms
,
589 fr
->x
, fr
->bV
? fr
->v
: nullptr, fr
->box
);
593 write_pdbfile(gmx_fio_getfp(status
->fio
), title
,
594 fr
->atoms
, fr
->x
, fr
->bPBC
? fr
->ePBC
: -1, fr
->box
,
599 write_g96_conf(gmx_fio_getfp(status
->fio
), title
, fr
, -1, nullptr);
602 gmx_fatal(FARGS
, "Sorry, write_trxframe can not write %s",
603 ftp2ext(gmx_fio_getftp(status
->fio
)));
609 int write_trx(t_trxstatus
*status
, int nind
, const int *ind
, const t_atoms
*atoms
,
610 int step
, real time
, matrix box
, rvec x
[], rvec
*v
,
615 clear_trxframe(&fr
, TRUE
);
620 fr
.bAtoms
= atoms
!= nullptr;
621 fr
.atoms
= const_cast<t_atoms
*>(atoms
);
624 fr
.bV
= v
!= nullptr;
627 copy_mat(box
, fr
.box
);
629 return write_trxframe_indexed(status
, &fr
, nind
, ind
, gc
);
632 void close_trx(t_trxstatus
*status
)
634 if (status
== nullptr)
638 gmx_tng_close(&status
->tng
);
641 gmx_fio_close(status
->fio
);
643 sfree(status
->persistent_line
);
645 sfree(status
->vmdplugin
);
647 /* The memory in status->xframe is lost here,
648 * but the read_first_x/read_next_x functions are deprecated anyhow.
649 * read_first_frame/read_next_frame and close_trx should be used.
654 t_trxstatus
*open_trx(const char *outfile
, const char *filemode
)
657 if (filemode
[0] != 'w' && filemode
[0] != 'a' && filemode
[1] != '+')
659 gmx_fatal(FARGS
, "Sorry, write_trx can only write");
665 stat
->fio
= gmx_fio_open(outfile
, filemode
);
669 static gmx_bool
gmx_next_frame(t_trxstatus
*status
, t_trxframe
*fr
)
676 if (gmx_trr_read_frame_header(status
->fio
, &sh
, &bOK
))
678 fr
->bDouble
= sh
.bDouble
;
679 fr
->natoms
= sh
.natoms
;
685 fr
->bFepState
= TRUE
;
686 fr
->lambda
= sh
.lambda
;
687 fr
->bBox
= sh
.box_size
> 0;
688 if (status
->flags
& (TRX_READ_X
| TRX_NEED_X
))
690 if (fr
->x
== nullptr)
692 snew(fr
->x
, sh
.natoms
);
694 fr
->bX
= sh
.x_size
> 0;
696 if (status
->flags
& (TRX_READ_V
| TRX_NEED_V
))
698 if (fr
->v
== nullptr)
700 snew(fr
->v
, sh
.natoms
);
702 fr
->bV
= sh
.v_size
> 0;
704 if (status
->flags
& (TRX_READ_F
| TRX_NEED_F
))
706 if (fr
->f
== nullptr)
708 snew(fr
->f
, sh
.natoms
);
710 fr
->bF
= sh
.f_size
> 0;
712 if (gmx_trr_read_frame_data(status
->fio
, &sh
, fr
->box
, fr
->x
, fr
->v
, fr
->f
))
718 fr
->not_ok
= DATA_NOT_OK
;
723 fr
->not_ok
= HEADER_NOT_OK
;
729 static gmx_bool
pdb_next_x(t_trxstatus
*status
, FILE *fp
, t_trxframe
*fr
)
734 // Initiate model_nr to -1 rather than NOTSET.
735 // It is not worthwhile introducing extra variables in the
736 // read_pdbfile call to verify that a model_nr was read.
737 int ePBC
, model_nr
= -1, na
;
738 char title
[STRLEN
], *time
, *step
;
741 atoms
.nr
= fr
->natoms
;
742 atoms
.atom
= nullptr;
743 atoms
.pdbinfo
= nullptr;
744 /* the other pointers in atoms should not be accessed if these are NULL */
747 na
= read_pdbfile(fp
, title
, &model_nr
, &atoms
, symtab
, fr
->x
, &ePBC
, boxpdb
, TRUE
, nullptr);
750 set_trxframe_ePBC(fr
, ePBC
);
751 if (nframes_read(status
) == 0)
753 fprintf(stderr
, " '%s', %d atoms\n", title
, fr
->natoms
);
758 fr
->bBox
= (boxpdb
[XX
][XX
] != 0.0);
761 copy_mat(boxpdb
, fr
->box
);
765 step
= std::strstr(title
, " step= ");
766 fr
->bStep
= ((step
!= nullptr) && sscanf(step
+7, "%" SCNd64
, &fr
->step
) == 1);
769 time
= std::strstr(title
, " t= ");
770 fr
->bTime
= ((time
!= nullptr) && sscanf(time
+4, "%lf", &dbl
) == 1);
779 if (na
!= fr
->natoms
)
781 gmx_fatal(FARGS
, "Number of atoms in pdb frame %d is %d instead of %d",
782 nframes_read(status
), na
, fr
->natoms
);
788 static int pdb_first_x(t_trxstatus
*status
, FILE *fp
, t_trxframe
*fr
)
792 fprintf(stderr
, "Reading frames from pdb file");
794 get_pdb_coordnum(fp
, &fr
->natoms
);
797 gmx_fatal(FARGS
, "\nNo coordinates in pdb file\n");
800 snew(fr
->x
, fr
->natoms
);
801 pdb_next_x(status
, fp
, fr
);
806 bool read_next_frame(const gmx_output_env_t
*oenv
, t_trxstatus
*status
, t_trxframe
*fr
)
810 gmx_bool bOK
, bMissingData
= FALSE
, bSkip
= FALSE
;
818 clear_trxframe(fr
, FALSE
);
822 /* Special treatment for TNG files */
827 ftp
= gmx_fio_getftp(status
->fio
);
832 bRet
= gmx_next_frame(status
, fr
);
835 /* Checkpoint files can not contain mulitple frames */
839 t_symtab
*symtab
= nullptr;
840 read_g96_conf(gmx_fio_getfp(status
->fio
), nullptr, nullptr, fr
,
841 symtab
, status
->persistent_line
);
842 bRet
= (fr
->natoms
> 0);
846 if (bTimeSet(TBEGIN
) && (status
->tf
< rTimeValue(TBEGIN
)))
848 if (xtc_seek_time(status
->fio
, rTimeValue(TBEGIN
), fr
->natoms
, TRUE
))
850 gmx_fatal(FARGS
, "Specified frame (time %f) doesn't exist or file corrupt/inconsistent.",
855 bRet
= (read_next_xtc(status
->fio
, fr
->natoms
, &fr
->step
, &fr
->time
, fr
->box
,
856 fr
->x
, &fr
->prec
, &bOK
) != 0);
857 fr
->bPrec
= (bRet
&& fr
->prec
> 0);
864 /* Actually the header could also be not ok,
865 but from bOK from read_next_xtc this can't be distinguished */
866 fr
->not_ok
= DATA_NOT_OK
;
870 bRet
= gmx_read_next_tng_frame(status
->tng
, fr
, nullptr, 0);
873 bRet
= pdb_next_x(status
, gmx_fio_getfp(status
->fio
), fr
);
876 bRet
= gro_next_x_or_v(gmx_fio_getfp(status
->fio
), fr
);
880 bRet
= read_next_vmd_frame(status
->vmdplugin
, fr
);
882 gmx_fatal(FARGS
, "DEATH HORROR in read_next_frame ftp=%s,status=%s",
883 ftp2ext(gmx_fio_getftp(status
->fio
)),
884 gmx_fio_getname(status
->fio
));
887 status
->tf
= fr
->time
;
891 bMissingData
= ((((status
->flags
& TRX_NEED_X
) != 0) && !fr
->bX
) ||
892 (((status
->flags
& TRX_NEED_V
) != 0) && !fr
->bV
) ||
893 (((status
->flags
& TRX_NEED_F
) != 0) && !fr
->bF
));
897 ct
= check_times2(fr
->time
, status
->t0
, fr
->bDouble
);
898 if (ct
== 0 || ((status
->flags
& TRX_DONT_SKIP
) && ct
< 0))
900 printcount(status
, oenv
, fr
->time
, FALSE
);
908 printcount(status
, oenv
, fr
->time
, TRUE
);
915 while (bRet
&& (bMissingData
|| bSkip
));
919 printlast(status
, oenv
, pt
);
922 printincomp(status
, fr
);
929 bool read_first_frame(const gmx_output_env_t
*oenv
, t_trxstatus
**status
,
930 const char *fn
, t_trxframe
*fr
, int flags
)
932 t_fileio
*fio
= nullptr;
933 gmx_bool bFirst
, bOK
;
934 int ftp
= fn2ftp(fn
);
936 clear_trxframe(fr
, TRUE
);
942 status_init( *status
);
944 (*status
)->flags
= flags
;
948 /* Special treatment for TNG files */
949 gmx_tng_open(fn
, 'r', &(*status
)->tng
);
953 fio
= (*status
)->fio
= gmx_fio_open(fn
, "r");
960 read_checkpoint_trxframe(fio
, fr
);
965 /* Can not rewind a compressed file, so open it twice */
966 if (!(*status
)->persistent_line
)
968 /* allocate the persistent line */
969 snew((*status
)->persistent_line
, STRLEN
+1);
971 t_symtab
*symtab
= nullptr;
972 read_g96_conf(gmx_fio_getfp(fio
), fn
, nullptr, fr
, symtab
, (*status
)->persistent_line
);
974 clear_trxframe(fr
, FALSE
);
975 if (flags
& (TRX_READ_X
| TRX_NEED_X
))
977 snew(fr
->x
, fr
->natoms
);
979 if (flags
& (TRX_READ_V
| TRX_NEED_V
))
981 snew(fr
->v
, fr
->natoms
);
983 (*status
)->fio
= gmx_fio_open(fn
, "r");
987 if (read_first_xtc(fio
, &fr
->natoms
, &fr
->step
, &fr
->time
, fr
->box
, &fr
->x
,
988 &fr
->prec
, &bOK
) == 0)
990 GMX_RELEASE_ASSERT(!bOK
, "Inconsistent results - OK status from read_first_xtc, but 0 atom coords read");
991 fr
->not_ok
= DATA_NOT_OK
;
996 printincomp(*status
, fr
);
1000 fr
->bPrec
= (fr
->prec
> 0);
1005 printcount(*status
, oenv
, fr
->time
, FALSE
);
1011 if (!gmx_read_next_tng_frame((*status
)->tng
, fr
, nullptr, 0))
1013 fr
->not_ok
= DATA_NOT_OK
;
1015 printincomp(*status
, fr
);
1019 printcount(*status
, oenv
, fr
->time
, FALSE
);
1024 pdb_first_x(*status
, gmx_fio_getfp(fio
), fr
);
1027 printcount(*status
, oenv
, fr
->time
, FALSE
);
1032 if (gro_first_x_or_v(gmx_fio_getfp(fio
), fr
))
1034 printcount(*status
, oenv
, fr
->time
, FALSE
);
1040 fprintf(stderr
, "The file format of %s is not a known trajectory format to GROMACS.\n"
1041 "Please make sure that the file is a trajectory!\n"
1042 "GROMACS will now assume it to be a trajectory and will try to open it using the VMD plug-ins.\n"
1043 "This will only work in case the VMD plugins are found and it is a trajectory format supported by VMD.\n", fn
);
1044 gmx_fio_fp_close(fio
); /*only close the file without removing FIO entry*/
1045 if (!read_first_vmd_frame(fn
, &(*status
)->vmdplugin
, fr
))
1047 gmx_fatal(FARGS
, "Not supported in read_first_frame: %s", fn
);
1050 gmx_fatal(FARGS
, "Not supported in read_first_frame: %s. Please make sure that the file is a trajectory.\n"
1051 "GROMACS is not compiled with plug-in support. Thus it cannot read non-GROMACS trajectory formats using the VMD plug-ins.\n"
1052 "Please compile with plug-in support if you want to read non-GROMACS trajectory formats.\n", fn
);
1055 (*status
)->tf
= fr
->time
;
1057 /* Return FALSE if we read a frame that's past the set ending time. */
1058 if (!bFirst
&& (!(flags
& TRX_DONT_SKIP
) && check_times(fr
->time
) > 0))
1060 (*status
)->t0
= fr
->time
;
1065 (!(flags
& TRX_DONT_SKIP
) && check_times(fr
->time
) < 0))
1067 /* Read a frame when no frame was read or the first was skipped */
1068 if (!read_next_frame(oenv
, *status
, fr
))
1073 (*status
)->t0
= fr
->time
;
1075 /* We need the number of atoms for random-access XTC searching, even when
1076 * we don't have access to the actual frame data.
1078 (*status
)->natoms
= fr
->natoms
;
1080 return (fr
->natoms
> 0);
1083 /***** C O O R D I N A T E S T U F F *****/
1085 int read_first_x(const gmx_output_env_t
*oenv
, t_trxstatus
**status
, const char *fn
,
1086 real
*t
, rvec
**x
, matrix box
)
1090 read_first_frame(oenv
, status
, fn
, &fr
, TRX_NEED_X
);
1092 snew((*status
)->xframe
, 1);
1093 (*(*status
)->xframe
) = fr
;
1094 *t
= (*status
)->xframe
->time
;
1095 *x
= (*status
)->xframe
->x
;
1096 copy_mat((*status
)->xframe
->box
, box
);
1098 return (*status
)->xframe
->natoms
;
1101 gmx_bool
read_next_x(const gmx_output_env_t
*oenv
, t_trxstatus
*status
, real
*t
,
1102 rvec x
[], matrix box
)
1106 status
->xframe
->x
= x
;
1107 /*xframe[status].x = x;*/
1108 bRet
= read_next_frame(oenv
, status
, status
->xframe
);
1109 *t
= status
->xframe
->time
;
1110 copy_mat(status
->xframe
->box
, box
);
1115 void rewind_trj(t_trxstatus
*status
)
1119 gmx_fio_rewind(status
->fio
);
1122 /***** T O P O L O G Y S T U F F ******/
1124 t_topology
*read_top(const char *fn
, int *ePBC
)
1130 epbc
= read_tpx_top(fn
, nullptr, nullptr, &natoms
, nullptr, nullptr, top
);