mdtools/md_code/InitializeFromFile.cpp

#include <iostream>
#include <cmath>

#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <unistd.h>
#include <sys/types.h>
#include <sys/stat.h>

#include "ReadWrite.hpp"
#include "simError.h"


InitializeFromFile :: InitializeFromFile( char *in_name ){

  c_in_file = fopen(in_name, "r");
  if(c_in_file == NULL){
    sprintf(painCave.errMsg,
            "Cannot open file: %s\n", in_name);
    painCave.isFatal = 1;
    simError();
  }
  
  strcpy( c_in_name, in_name);
  return;  
}

InitializeFromFile :: ~InitializeFromFile( ){
  
  int error;
  error = fclose( c_in_file );
  if( error ){
    sprintf( painCave.errMsg,
             "Error closing %s\n", c_in_name );
    simError();
  }
  return;
}


void InitializeFromFile :: read_xyz( SimInfo* entry_plug ){

  int i; // loop counter

  int n_atoms; // the number of atoms
  char read_buffer[2000]; //the line buffer for reading 
  char *eof_test; // ptr to see when we reach the end of the file 
  char *foo; // the pointer to the current string token 

  double rx, ry, rz; // position place holders
  double vx, vy, vz; // velocity placeholders
  double q[4]; // the quaternions
  double jx, jy, jz; // angular velocity placeholders;
  double qSqr, qLength; // needed to normalize the quaternion vector.

  eof_test = fgets(read_buffer, sizeof(read_buffer), c_in_file);
  if( eof_test == NULL ){
    std::cerr << "error reading 1st line of" << c_in_name << "\n";
  }

  (void)sscanf(read_buffer, "%d", &n_atoms);

  Atom **atoms = entry_plug->atoms;
  DirectionalAtom* dAtom;

  if( n_atoms != entry_plug->n_atoms ){
    sprintf( painCave.errMsg,
             "Initialize from File error. %s n_atoms, %d, "
             "does not match the BASS file's n_atoms, %d.\n",
             c_in_name, n_atoms, entry_plug->n_atoms );
    painCave.isFatal = 1;
    simError();
  }
  
  //read and toss the comment line 
  
  eof_test = fgets(read_buffer, sizeof(read_buffer), c_in_file);
  if(eof_test == NULL){
    sprintf( painCave.errMsg,
             "error in reading commment in %s\n", c_in_name);
    painCave.isFatal = 1;
    simError();
  }

  for( i=0; i < n_atoms; i++){
    
    eof_test = fgets(read_buffer, sizeof(read_buffer), c_in_file);
    if(eof_test == NULL){
      sprintf(painCave.errMsg,
              "error in reading file %s\n"
              "natoms  = %d; index = %d\n"
              "error reading the line from the file.\n",
              c_in_name, n_atoms, i );
      painCave.isFatal = 1;
      simError();
    }

    foo = strtok(read_buffer, " ,;\t");
    
    // check the atom name to the current atom

    if( strcmp( foo, atoms[i]->getType() ) ){
      sprintf( painCave.errMsg,
               "Initialize from file error. Atom %s at index %d "
               "in file %s does not"
               " match the BASS atom %s.\n",
               foo, i, c_in_name, atoms[i]->getType() );
      painCave.isFatal = 1;
      simError();
    }
    
    // get the positions

    foo = strtok(NULL, " ,;\t");
    if(foo == NULL){
      sprintf( painCave.errMsg,
               "error in reading postition x from %s\n"
               "natoms  = %d, index = %d\n",
               c_in_name, n_atoms, i );
      painCave.isFatal = 1;
      simError();
    }
    (void)sscanf( foo, "%lf", &rx );
   
    foo = strtok(NULL, " ,;\t");
    if(foo == NULL){
      sprintf( painCave.errMsg,
               "error in reading postition y from %s\n"
               "natoms  = %d, index = %d\n",
               c_in_name, n_atoms, i );
      painCave.isFatal = 1;
      simError();
    }
    (void)sscanf( foo, "%lf", &ry );
     
    foo = strtok(NULL, " ,;\t");
    if(foo == NULL){
      sprintf( painCave.errMsg,
               "error in reading postition z from %s\n"
               "natoms  = %d, index = %d\n",
               c_in_name, n_atoms, i );
      painCave.isFatal = 1;
      simError();
    }
    (void)sscanf( foo, "%lf", &rz );
    
    // get the velocities

    foo = strtok(NULL, " ,;\t");
    if(foo == NULL){
      sprintf( painCave.errMsg,
               "error in reading velocity x from %s\n"
               "natoms  = %d, index = %d\n",
               c_in_name, n_atoms, i );
      painCave.isFatal = 1;
      simError();
    }
    (void)sscanf( foo, "%lf", &vx );
    
    foo = strtok(NULL, " ,;\t");
    if(foo == NULL){
      sprintf( painCave.errMsg,
               "error in reading velocity y from %s\n"
               "natoms  = %d, index = %d\n",
               c_in_name, n_atoms, i );
      painCave.isFatal = 1;
      simError();
    }
    (void)sscanf( foo, "%lf", &vy );
    
    foo = strtok(NULL, " ,;\t");
    if(foo == NULL){
      sprintf( painCave.errMsg,
               "error in reading velocity z from %s\n"
               "natoms  = %d, index = %d\n",
               c_in_name, n_atoms, i );
      painCave.isFatal = 1;
      simError();
    }
    (void)sscanf( foo, "%lf", &vz );
    
    
    // get the quaternions
    
    if( atoms[i]->isDirectional() ){
      
      foo = strtok(NULL, " ,;\t");
      if(foo == NULL){
        sprintf(painCave.errMsg,
                "error in reading quaternion 0 from %s\n"
                "natoms  = %d, index = %d\n",
                c_in_name, n_atoms, i );
        painCave.isFatal = 1;
        simError();
      }
      (void)sscanf( foo, "%lf", &q[0] );
      
      foo = strtok(NULL, " ,;\t");
      if(foo == NULL){
        sprintf( painCave.errMsg,
                 "error in reading quaternion 1 from %s\n"
                 "natoms  = %d, index = %d\n",
                 c_in_name, n_atoms, i );
        painCave.isFatal = 1;
        simError();
      }
      (void)sscanf( foo, "%lf", &q[1] );
      
      foo = strtok(NULL, " ,;\t");
      if(foo == NULL){
        sprintf( painCave.errMsg,
                 "error in reading quaternion 2 from %s\n"
                 "natoms  = %d, index = %d\n",
                 c_in_name, n_atoms, i );
        painCave.isFatal = 1;
        simError();
      }
      (void)sscanf( foo, "%lf", &q[2] );
      
      foo = strtok(NULL, " ,;\t");
      if(foo == NULL){
        sprintf( painCave.errMsg,
                 "error in reading quaternion 3 from %s\n"
                 "natoms  = %d, index = %d\n",
                 c_in_name, n_atoms, i );
        painCave.isFatal = 1;
        simError();
      }
      (void)sscanf( foo, "%lf", &q[3] );
      
      // get the angular velocities
      
      foo = strtok(NULL, " ,;\t");
      if(foo == NULL){
        sprintf( painCave.errMsg,
                 "error in reading angular momentum jx from %s\n"
                 "natoms  = %d, index = %d\n",
                 c_in_name, n_atoms, i );
        painCave.isFatal = 1;
        simError();
      }
      (void)sscanf( foo, "%lf", &jx );
      
      foo = strtok(NULL, " ,;\t");
      if(foo == NULL){
        sprintf( painCave.errMsg,
                 "error in reading angular momentum jy from %s\n"
                 "natoms  = %d, index = %d\n",
                 c_in_name, n_atoms, i );
        painCave.isFatal = 1;
        simError();
      }
      (void)sscanf( foo, "%lf", &jy );
      
      foo = strtok(NULL, " ,;\t");
      if(foo == NULL){
        sprintf( painCave.errMsg,
                 "error in reading angular momentum jz from %s\n"
                 "natoms  = %d, index = %d\n",
                 c_in_name, n_atoms, i );
        painCave.isFatal = 1;
        simError();
      }
      (void)sscanf( foo, "%lf", &jz );
      
      dAtom = ( DirectionalAtom* )atoms[i];

      // check that the quaternion vector is normalized

      qSqr = (q[0] * q[0]) + (q[1] * q[1]) + (q[2] * q[2]) + (q[3] * q[3]);
      
      qLength = sqrt( qSqr );
      q[0] = q[0] / qLength;
      q[1] = q[1] / qLength;
      q[2] = q[2] / qLength;
      q[3] = q[3] / qLength;
      
      dAtom->setQ( q );
      
      // add the angular velocities

      dAtom->setJx( jx );
      dAtom->setJy( jy );
      dAtom->setJz( jz );
    }
    
    // add the positions and velocities to the atom
    
    atoms[i]->setX( rx );
    atoms[i]->setY( ry );
    atoms[i]->setZ( rz );
    
    atoms[i]->set_vx( vx );
    atoms[i]->set_vy( vy );
    atoms[i]->set_vz( vz );
    
  }
}
Revision:	162
Committed:	Thu Oct 31 21:20:49 2002 UTC (21 years, 8 months ago) by mmeineke
File size:	7546 byte(s)
Log Message:	* empty log message *
#	User	Rev	Content
1	mmeineke	10	#include <iostream>
2			#include <cmath>
3
4			#include <stdio.h>
5			#include <stdlib.h>
6			#include <string.h>
7			#include <unistd.h>
8			#include <sys/types.h>
9			#include <sys/stat.h>
10
11			#include "ReadWrite.hpp"
12	mmeineke	162	#include "simError.h"
13	mmeineke	10
14
15			InitializeFromFile :: InitializeFromFile( char *in_name ){
16
17			c_in_file = fopen(in_name, "r");
18			if(c_in_file == NULL){
19	mmeineke	162	sprintf(painCave.errMsg,
20			"Cannot open file: %s\n", in_name);
21			painCave.isFatal = 1;
22			simError();
23	mmeineke	10	}
24
25			strcpy( c_in_name, in_name);
26			return;
27			}
28
29			InitializeFromFile :: ~InitializeFromFile( ){
30
31			int error;
32			error = fclose( c_in_file );
33			if( error ){
34	mmeineke	162	sprintf( painCave.errMsg,
35			"Error closing %s\n", c_in_name );
36			simError();
37	mmeineke	10	}
38			return;
39			}
40
41
42			void InitializeFromFile :: read_xyz( SimInfo* entry_plug ){
43
44			int i; // loop counter
45
46			int n_atoms; // the number of atoms
47			char read_buffer[2000]; //the line buffer for reading
48			char *eof_test; // ptr to see when we reach the end of the file
49			char *foo; // the pointer to the current string token
50
51			double rx, ry, rz; // position place holders
52			double vx, vy, vz; // velocity placeholders
53			double q[4]; // the quaternions
54			double jx, jy, jz; // angular velocity placeholders;
55			double qSqr, qLength; // needed to normalize the quaternion vector.
56
57			eof_test = fgets(read_buffer, sizeof(read_buffer), c_in_file);
58			if( eof_test == NULL ){
59			std::cerr << "error reading 1st line of" << c_in_name << "\n";
60			}
61
62			(void)sscanf(read_buffer, "%d", &n_atoms);
63
64			Atom **atoms = entry_plug->atoms;
65			DirectionalAtom* dAtom;
66
67			if( n_atoms != entry_plug->n_atoms ){
68	mmeineke	162	sprintf( painCave.errMsg,
69	mmeineke	10	"Initialize from File error. %s n_atoms, %d, "
70			"does not match the BASS file's n_atoms, %d.\n",
71			c_in_name, n_atoms, entry_plug->n_atoms );
72	mmeineke	162	painCave.isFatal = 1;
73			simError();
74	mmeineke	10	}
75
76			//read and toss the comment line
77
78			eof_test = fgets(read_buffer, sizeof(read_buffer), c_in_file);
79			if(eof_test == NULL){
80	mmeineke	162	sprintf( painCave.errMsg,
81			"error in reading commment in %s\n", c_in_name);
82			painCave.isFatal = 1;
83			simError();
84	mmeineke	10	}
85
86			for( i=0; i < n_atoms; i++){
87
88			eof_test = fgets(read_buffer, sizeof(read_buffer), c_in_file);
89			if(eof_test == NULL){
90	mmeineke	162	sprintf(painCave.errMsg,
91			"error in reading file %s\n"
92			"natoms = %d; index = %d\n"
93			"error reading the line from the file.\n",
94			c_in_name, n_atoms, i );
95			painCave.isFatal = 1;
96			simError();
97	mmeineke	10	}
98
99			foo = strtok(read_buffer, " ,;\t");
100
101			// check the atom name to the current atom
102
103			if( strcmp( foo, atoms[i]->getType() ) ){
104	mmeineke	162	sprintf( painCave.errMsg,
105	mmeineke	10	"Initialize from file error. Atom %s at index %d "
106			"in file %s does not"
107			" match the BASS atom %s.\n",
108			foo, i, c_in_name, atoms[i]->getType() );
109	mmeineke	162	painCave.isFatal = 1;
110			simError();
111	mmeineke	10	}
112
113			// get the positions
114
115			foo = strtok(NULL, " ,;\t");
116			if(foo == NULL){
117	mmeineke	162	sprintf( painCave.errMsg,
118			"error in reading postition x from %s\n"
119			"natoms = %d, index = %d\n",
120			c_in_name, n_atoms, i );
121			painCave.isFatal = 1;
122			simError();
123	mmeineke	10	}
124			(void)sscanf( foo, "%lf", &rx );
125
126			foo = strtok(NULL, " ,;\t");
127			if(foo == NULL){
128	mmeineke	162	sprintf( painCave.errMsg,
129			"error in reading postition y from %s\n"
130			"natoms = %d, index = %d\n",
131			c_in_name, n_atoms, i );
132			painCave.isFatal = 1;
133			simError();
134	mmeineke	10	}
135			(void)sscanf( foo, "%lf", &ry );
136
137			foo = strtok(NULL, " ,;\t");
138			if(foo == NULL){
139	mmeineke	162	sprintf( painCave.errMsg,
140			"error in reading postition z from %s\n"
141			"natoms = %d, index = %d\n",
142			c_in_name, n_atoms, i );
143			painCave.isFatal = 1;
144			simError();
145	mmeineke	10	}
146			(void)sscanf( foo, "%lf", &rz );
147
148			// get the velocities
149
150			foo = strtok(NULL, " ,;\t");
151			if(foo == NULL){
152	mmeineke	162	sprintf( painCave.errMsg,
153			"error in reading velocity x from %s\n"
154			"natoms = %d, index = %d\n",
155			c_in_name, n_atoms, i );
156			painCave.isFatal = 1;
157			simError();
158	mmeineke	10	}
159			(void)sscanf( foo, "%lf", &vx );
160
161			foo = strtok(NULL, " ,;\t");
162			if(foo == NULL){
163	mmeineke	162	sprintf( painCave.errMsg,
164			"error in reading velocity y from %s\n"
165			"natoms = %d, index = %d\n",
166			c_in_name, n_atoms, i );
167			painCave.isFatal = 1;
168			simError();
169	mmeineke	10	}
170			(void)sscanf( foo, "%lf", &vy );
171
172			foo = strtok(NULL, " ,;\t");
173			if(foo == NULL){
174	mmeineke	162	sprintf( painCave.errMsg,
175			"error in reading velocity z from %s\n"
176			"natoms = %d, index = %d\n",
177			c_in_name, n_atoms, i );
178			painCave.isFatal = 1;
179			simError();
180	mmeineke	10	}
181			(void)sscanf( foo, "%lf", &vz );
182
183
184			// get the quaternions
185
186			if( atoms[i]->isDirectional() ){
187
188			foo = strtok(NULL, " ,;\t");
189			if(foo == NULL){
190	mmeineke	162	sprintf(painCave.errMsg,
191			"error in reading quaternion 0 from %s\n"
192			"natoms = %d, index = %d\n",
193			c_in_name, n_atoms, i );
194			painCave.isFatal = 1;
195			simError();
196	mmeineke	10	}
197			(void)sscanf( foo, "%lf", &q[0] );
198
199			foo = strtok(NULL, " ,;\t");
200			if(foo == NULL){
201	mmeineke	162	sprintf( painCave.errMsg,
202			"error in reading quaternion 1 from %s\n"
203			"natoms = %d, index = %d\n",
204			c_in_name, n_atoms, i );
205			painCave.isFatal = 1;
206			simError();
207	mmeineke	10	}
208			(void)sscanf( foo, "%lf", &q[1] );
209
210			foo = strtok(NULL, " ,;\t");
211			if(foo == NULL){
212	mmeineke	162	sprintf( painCave.errMsg,
213			"error in reading quaternion 2 from %s\n"
214			"natoms = %d, index = %d\n",
215			c_in_name, n_atoms, i );
216			painCave.isFatal = 1;
217			simError();
218	mmeineke	10	}
219			(void)sscanf( foo, "%lf", &q[2] );
220	mmeineke	162
221	mmeineke	10	foo = strtok(NULL, " ,;\t");
222			if(foo == NULL){
223	mmeineke	162	sprintf( painCave.errMsg,
224			"error in reading quaternion 3 from %s\n"
225			"natoms = %d, index = %d\n",
226			c_in_name, n_atoms, i );
227			painCave.isFatal = 1;
228			simError();
229	mmeineke	10	}
230			(void)sscanf( foo, "%lf", &q[3] );
231
232			// get the angular velocities
233
234			foo = strtok(NULL, " ,;\t");
235			if(foo == NULL){
236	mmeineke	162	sprintf( painCave.errMsg,
237			"error in reading angular momentum jx from %s\n"
238			"natoms = %d, index = %d\n",
239			c_in_name, n_atoms, i );
240			painCave.isFatal = 1;
241			simError();
242	mmeineke	10	}
243			(void)sscanf( foo, "%lf", &jx );
244
245			foo = strtok(NULL, " ,;\t");
246			if(foo == NULL){
247	mmeineke	162	sprintf( painCave.errMsg,
248			"error in reading angular momentum jy from %s\n"
249			"natoms = %d, index = %d\n",
250			c_in_name, n_atoms, i );
251			painCave.isFatal = 1;
252			simError();
253	mmeineke	10	}
254			(void)sscanf( foo, "%lf", &jy );
255
256			foo = strtok(NULL, " ,;\t");
257			if(foo == NULL){
258	mmeineke	162	sprintf( painCave.errMsg,
259			"error in reading angular momentum jz from %s\n"
260			"natoms = %d, index = %d\n",
261			c_in_name, n_atoms, i );
262			painCave.isFatal = 1;
263			simError();
264	mmeineke	10	}
265			(void)sscanf( foo, "%lf", &jz );
266
267			dAtom = ( DirectionalAtom* )atoms[i];
268
269			// check that the quaternion vector is normalized
270
271			qSqr = (q[0] * q[0]) + (q[1] * q[1]) + (q[2] * q[2]) + (q[3] * q[3]);
272
273			qLength = sqrt( qSqr );
274			q[0] = q[0] / qLength;
275			q[1] = q[1] / qLength;
276			q[2] = q[2] / qLength;
277			q[3] = q[3] / qLength;
278
279			dAtom->setQ( q );
280
281			// add the angular velocities
282
283			dAtom->setJx( jx );
284			dAtom->setJy( jy );
285			dAtom->setJz( jz );
286			}
287
288			// add the positions and velocities to the atom
289
290			atoms[i]->setX( rx );
291			atoms[i]->setY( ry );
292			atoms[i]->setZ( rz );
293
294			atoms[i]->set_vx( vx );
295			atoms[i]->set_vy( vy );
296			atoms[i]->set_vz( vz );
297
298			}
299			}