src/gmxlib/xtcio.c

   1 /*
   2  * This file is part of the GROMACS molecular simulation package.
   3  *
   4  * Copyright (c) 1991-2000, University of Groningen, The Netherlands.
   5  * Copyright (c) 2001-2004, The GROMACS development team,
   6  * check out http://www.gromacs.org for more information.
   7  * Copyright (c) 2012, by the GROMACS development team, led by
   8  * David van der Spoel, Berk Hess, Erik Lindahl, and including many
   9  * others, as listed in the AUTHORS file in the top-level source
  10  * directory and at http://www.gromacs.org.
  11  *
  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.
  16  *
  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.
  21  *
  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.
  26  *
  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.
  34  *
  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.
  37  */
  38 #ifdef HAVE_CONFIG_H
  39 #include <config.h>
  40 #endif
  41
  42 #include <string.h>
  43 #include "typedefs.h"
  44 #include "xdrf.h"
  45 #include "gmxfio.h"
  46 #include "xtcio.h"
  47 #include "smalloc.h"
  48 #include "vec.h"
  49 #include "futil.h"
  50 #include "gmx_fatal.h"
  51
  52 #define XTC_MAGIC 1995
  53
  54
  55 static int xdr_r2f(XDR *xdrs,real *r,gmx_bool bRead)
  56 {
  57 #ifdef GMX_DOUBLE
  58     float f;
  59     int   ret;
  60
  61     if (!bRead)
  62       f = *r;
  63     ret = xdr_float(xdrs,&f);
  64     if (bRead)
  65       *r = f;
  66
  67     return ret;
  68 #else
  69     return xdr_float(xdrs,(float *)r);
  70 #endif
  71 }
  72
  73
  74 t_fileio *open_xtc(const char *fn,const char *mode)
  75 {
  76     return gmx_fio_open(fn,mode);
  77 }
  78
  79 void close_xtc(t_fileio *fio)
  80 {
  81     gmx_fio_close(fio);
  82 }
  83
  84 static void check_xtc_magic(int magic)
  85 {
  86   if (magic != XTC_MAGIC)
  87     gmx_fatal(FARGS,"Magic Number Error in XTC file (read %d, should be %d)",
  88                 magic,XTC_MAGIC);
  89 }
  90
  91 int xtc_check(const char *str,gmx_bool bResult,const char *file,int line)
  92 {
  93   if (!bResult) {
  94     if (debug)
  95       fprintf(debug,"\nXTC error: read/write of %s failed, "
  96               "source file %s, line %d\n",str,file,line);
  97     return 0;
  98   }
  99   return 1;
 100 }
 101
 102 void xtc_check_fat_err(const char *str,gmx_bool bResult,const char *file,int line)
 103 {
 104   if (!bResult) {
 105     gmx_fatal(FARGS,"XTC read/write of %s failed, "
 106                 "source file %s, line %d\n",str,file,line);
 107   }
 108 }
 109
 110 static int xtc_header(XDR *xd,int *magic,int *natoms,int *step,real *time,
 111                       gmx_bool bRead,gmx_bool *bOK)
 112 {
 113   int result;
 114
 115   if (xdr_int(xd,magic) == 0)
 116     return 0;
 117   result=XTC_CHECK("natoms", xdr_int(xd,natoms));  /* number of atoms */
 118   if (result)
 119     result=XTC_CHECK("step",   xdr_int(xd,step));    /* frame number    */
 120   if (result)
 121     result=XTC_CHECK("time",   xdr_r2f(xd,time,bRead));   /* time */
 122   *bOK=(result!=0);
 123
 124   return result;
 125 }
 126
 127 static int xtc_coord(XDR *xd,int *natoms,matrix box,rvec *x,real *prec, gmx_bool bRead)
 128 {
 129   int i,j,result;
 130 #ifdef GMX_DOUBLE
 131   float *ftmp;
 132   float fprec;
 133 #endif
 134
 135   /* box */
 136   result=1;
 137   for(i=0; ((i<DIM) && result); i++)
 138     for(j=0; ((j<DIM) && result); j++)
 139       result=XTC_CHECK("box",xdr_r2f(xd,&(box[i][j]),bRead));
 140
 141   if (!result)
 142       return result;
 143
 144 #ifdef GMX_DOUBLE
 145   /* allocate temp. single-precision array */
 146   snew(ftmp,(*natoms)*DIM);
 147
 148   /* Copy data to temp. array if writing */
 149   if(!bRead)
 150   {
 151       for(i=0; (i<*natoms); i++)
 152       {
 153           ftmp[DIM*i+XX]=x[i][XX];
 154           ftmp[DIM*i+YY]=x[i][YY];
 155           ftmp[DIM*i+ZZ]=x[i][ZZ];
 156       }
 157       fprec = *prec;
 158   }
 159   result=XTC_CHECK("x",xdr3dfcoord(xd,ftmp,natoms,&fprec));
 160
 161   /* Copy from temp. array if reading */
 162   if(bRead)
 163   {
 164       for(i=0; (i<*natoms); i++)
 165       {
 166           x[i][XX] = ftmp[DIM*i+XX];
 167           x[i][YY] = ftmp[DIM*i+YY];
 168           x[i][ZZ] = ftmp[DIM*i+ZZ];
 169       }
 170       *prec = fprec;
 171   }
 172   sfree(ftmp);
 173 #else
 174     result=XTC_CHECK("x",xdr3dfcoord(xd,x[0],natoms,prec));
 175 #endif
 176
 177   return result;
 178 }
 179
 180
 181
 182 int write_xtc(t_fileio *fio,
 183               int natoms,int step,real time,
 184               matrix box,rvec *x,real prec)
 185 {
 186   int magic_number = XTC_MAGIC;
 187   XDR *xd;
 188   gmx_bool bDum;
 189   int bOK;
 190
 191   xd = gmx_fio_getxdr(fio);
 192   /* write magic number and xtc identidier */
 193   if (xtc_header(xd,&magic_number,&natoms,&step,&time,FALSE,&bDum) == 0)
 194   {
 195           return 0;
 196   }
 197
 198   /* write data */
 199   bOK = xtc_coord(xd,&natoms,box,x,&prec,FALSE); /* bOK will be 1 if writing went well */
 200
 201   if(bOK)
 202   {
 203           if(gmx_fio_flush(fio) !=0)
 204           {
 205                   bOK = 0;
 206           }
 207   }
 208   return bOK;  /* 0 if bad, 1 if writing went well */
 209 }
 210
 211 int read_first_xtc(t_fileio *fio,int *natoms,int *step,real *time,
 212                    matrix box,rvec **x,real *prec,gmx_bool *bOK)
 213 {
 214   int magic;
 215   XDR *xd;
 216
 217   *bOK=TRUE;
 218   xd = gmx_fio_getxdr(fio);
 219
 220   /* read header and malloc x */
 221   if ( !xtc_header(xd,&magic,natoms,step,time,TRUE,bOK))
 222     return 0;
 223
 224   /* Check magic number */
 225   check_xtc_magic(magic);
 226
 227   snew(*x,*natoms);
 228
 229   *bOK=xtc_coord(xd,natoms,box,*x,prec,TRUE);
 230
 231   return *bOK;
 232 }
 233
 234 int read_next_xtc(t_fileio* fio,
 235                   int natoms,int *step,real *time,
 236                   matrix box,rvec *x,real *prec,gmx_bool *bOK)
 237 {
 238   int magic;
 239   int n;
 240   XDR *xd;
 241
 242   *bOK=TRUE;
 243   xd = gmx_fio_getxdr(fio);
 244
 245   /* read header */
 246   if (!xtc_header(xd,&magic,&n,step,time,TRUE,bOK))
 247     return 0;
 248
 249   /* Check magic number */
 250   check_xtc_magic(magic);
 251
 252   if (n > natoms) {
 253     gmx_fatal(FARGS, "Frame contains more atoms (%d) than expected (%d)",
 254               n, natoms);
 255   }
 256
 257   *bOK=xtc_coord(xd,&natoms,box,x,prec,TRUE);
 258
 259   return *bOK;
 260 }
 261
 262