ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/OpenMD/branches/development/src/io/RestReader.cpp
(Generate patch)

Comparing trunk/src/io/RestReader.cpp (file contents):
Revision 1407 by cli2, Wed Jan 20 16:04:40 2010 UTC vs.
Revision 1408 by gezelter, Fri Jan 22 21:26:29 2010 UTC

# Line 55 | Line 55
55   #include "primitives/Molecule.hpp"
56   #include "utils/simError.h"
57   #include "utils/MemoryUtils.hpp"
58 < #include "utils/StringTokenizer.hpp"
58 > #include "utils/StringTokenizer.hpp"the
59   #include "restraints/ObjectRestraint.hpp"
60   #include "restraints/MolecularRestraint.hpp"
61  
# Line 191 | Line 191 | namespace OpenMD {
191  
192      scanFile();
193  
194 –
194      readSet();
195  
196  
# Line 280 | Line 279 | namespace OpenMD {
279  
280      int index = tokenizer.nextTokenAsInt();
281  
282 +    StuntDouble* integrableObject = info_->getIOIndexToIntegrableObject(index);
283 +
284 +    if (integrableObject == NULL) {
285 +      return;
286 +    }
287    
288 <    if (index < 1116){
289 <      std::string type = tokenizer.nextToken();
290 <      
291 <      Vector3d pos;
292 <      
293 <      pos[0] = tokenizer.nextTokenAsDouble();
294 <      pos[1] = tokenizer.nextTokenAsDouble();
295 <      pos[2] = tokenizer.nextTokenAsDouble();
296 <      
297 <      all_pos_[index] = pos;
294 <    }else{
295 <      
296 <      bool exist = false;
297 <      int indexCount = 0;
298 <      
299 <      while ( (!exist) && (indexCount < stuntDoubleIndex_.size())){
300 <        if (index == stuntDoubleIndex_[indexCount])
301 <          exist = true;
302 <        indexCount++;
303 <      }
304 <      
305 <      StuntDouble* integrableObject;
306 <      
307 <      if (exist){
308 <        
309 <        integrableObject = info_->getIOIndexToIntegrableObject(index);
310 <        
311 <        int compareindex = integrableObject->getGlobalIntegrableObjectIndex();
312 <        
313 <        if (integrableObject == NULL) {
314 <          return;
315 <        }  
316 <        
317 <        std::string type = tokenizer.nextToken();
318 <        
319 <        Vector3d pos;
320 <        
288 >    std::string type = tokenizer.nextToken();
289 >    int size = type.size();
290 >    
291 >    Vector3d pos;
292 >    Quat4d q;
293 >
294 >    for(int i = 0; i < size; ++i) {
295 >      switch(type[i]) {
296 >
297 >      case 'p': {
298          pos[0] = tokenizer.nextTokenAsDouble();
299          pos[1] = tokenizer.nextTokenAsDouble();
300          pos[2] = tokenizer.nextTokenAsDouble();
301 <        
301 >        break;
302 >      }
303 >      case 'v' : {
304          Vector3d vel;
305          vel[0] = tokenizer.nextTokenAsDouble();
306          vel[1] = tokenizer.nextTokenAsDouble();
307          vel[2] = tokenizer.nextTokenAsDouble();
308 <        
309 <        
310 <        Quat4d q;
308 >        break;
309 >      }
310 >
311 >      case 'q' : {
312          if (integrableObject->isDirectional()) {
313            
314            q[0] = tokenizer.nextTokenAsDouble();
315            q[1] = tokenizer.nextTokenAsDouble();
316            q[2] = tokenizer.nextTokenAsDouble();
317            q[3] = tokenizer.nextTokenAsDouble();
318 <        }  
319 <        // keep the position in case we need it for a molecular restraint:
320 <        
321 <        all_pos_[index] = pos;
322 <        
323 <        // is this io restrained?
324 <        GenericData* data = integrableObject->getPropertyByName("Restraint");
325 <        ObjectRestraint* oRest;
318 >          
319 >          RealType qlen = q.length();
320 >          if (qlen < OpenMD::epsilon) { //check quaternion is not equal to 0
321 >            
322 >            sprintf(painCave.errMsg,
323 >                    "RestReader Error: initial quaternion error (q0^2 + q1^2 + q2^2 + q3^2) ~ 0\n");
324 >            painCave.isFatal = 1;
325 >            simError();            
326 >          }  
327 >              
328 >          q.normalize();
329 >        }          
330 >        break;
331 >      }  
332 >      case 'j' : {
333 >        Vector3d ji;
334 >        if (integrableObject->isDirectional()) {
335 >          ji[0] = tokenizer.nextTokenAsDouble();
336 >          ji[1] = tokenizer.nextTokenAsDouble();
337 >          ji[2] = tokenizer.nextTokenAsDouble();
338 >        }
339 >        break;
340 >      }  
341 >      case 'f': {        
342 >        Vector3d force;
343 >        force[0] = tokenizer.nextTokenAsDouble();
344 >        force[1] = tokenizer.nextTokenAsDouble();
345 >        force[2] = tokenizer.nextTokenAsDouble();          
346 >        break;
347 >      }
348 >      case 't' : {        
349 >        Vector3d torque;
350 >        torque[0] = tokenizer.nextTokenAsDouble();
351 >        torque[1] = tokenizer.nextTokenAsDouble();
352 >        torque[2] = tokenizer.nextTokenAsDouble();          
353 >        break;
354 >      }
355 >      default: {
356 >        sprintf(painCave.errMsg,
357 >                "RestReader Error: %s is an unrecognized type\n", type.c_str());
358 >        painCave.isFatal = 1;
359 >        simError();    
360 >        break;  
361 >      }
362 >      }
363 >      // keep the position in case we need it for a molecular restraint:
364 >
365 >      all_pos_[index] = pos;      
366          
367 <        if (data != NULL) {
368 <          // make sure we can reinterpret the generic data as restraint data:
369 <          RestraintData* restData= dynamic_cast<RestraintData*>(data);        
370 <          if (restData != NULL) {
371 <            // make sure we can reinterpet the restraint data as a pointer to
367 >      // is this io restrained?
368 >      GenericData* data = integrableObject->getPropertyByName("Restraint");
369 >      ObjectRestraint* oRest;
370 >      
371 >      if (data != NULL) {
372 >        // make sure we can reinterpret the generic data as restraint data:
373 >        RestraintData* restData= dynamic_cast<RestraintData*>(data);        
374 >        if (restData != NULL) {
375 >          // make sure we can reinterpet the restraint data as a pointer to
376              // an ObjectRestraint:
377 <            oRest = dynamic_cast<ObjectRestraint*>(restData->getData());
378 <            if (oRest != NULL) {          
379 <              if (integrableObject->isDirectional()) {
380 <                oRest->setReferenceStructure(pos, q.toRotationMatrix3());
381 <              } else {                          
382 <                oRest->setReferenceStructure(pos);
359 <              }
377 >          oRest = dynamic_cast<ObjectRestraint*>(restData->getData());
378 >          if (oRest != NULL) {          
379 >            if (integrableObject->isDirectional()) {
380 >              oRest->setReferenceStructure(pos, q.toRotationMatrix3());
381 >            } else {                          
382 >              oRest->setReferenceStructure(pos);
383              }
384            }
385          }

Diff Legend

– Removed lines
+ Added lines
< Changed lines
> Changed lines