ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/OpenMD/trunk/src/selection/SelectionEvaluator.cpp
(Generate patch)

Comparing trunk/src/selection/SelectionEvaluator.cpp (file contents):
Revision 1903 by gezelter, Tue Jul 16 18:58:08 2013 UTC vs.
Revision 2071 by gezelter, Sat Mar 7 21:41:51 2015 UTC

# Line 53 | Line 53 | namespace OpenMD {
53  
54   namespace OpenMD {
55  
56 –
56    SelectionEvaluator::SelectionEvaluator(SimInfo* si)
57      : info(si), nameFinder(info), distanceFinder(info), hullFinder(info),
58 <      indexFinder(info), hasSurfaceArea_(false),
59 <      isLoaded_(false){    
60 <      nStuntDouble = info->getNGlobalAtoms() + info->getNGlobalRigidBodies();
61 <    }            
62 <
58 >      indexFinder(info), isLoaded_(false), hasSurfaceArea_(false) {
59 >    nObjects.push_back(info->getNGlobalAtoms() + info->getNGlobalRigidBodies());
60 >    nObjects.push_back(info->getNGlobalBonds());
61 >    nObjects.push_back(info->getNGlobalBends());
62 >    nObjects.push_back(info->getNGlobalTorsions());
63 >    nObjects.push_back(info->getNGlobalInversions());
64 >    nObjects.push_back(info->getNGlobalMolecules());    
65 >  }
66 >  
67    bool SelectionEvaluator::loadScript(const std::string& filename,
68                                        const std::string& script) {
69      clearDefinitionsAndLoadPredefined();
# Line 127 | Line 130 | namespace OpenMD {
130      return loadScript(filename, script);
131    }
132    
133 <  void SelectionEvaluator::instructionDispatchLoop(OpenMDBitSet& bs){
133 >  void SelectionEvaluator::instructionDispatchLoop(SelectionSet& bs){
134      
135      while ( pc < aatoken.size()) {
136        statement = aatoken[pc++];
# Line 148 | Line 151 | namespace OpenMD {
151  
152    }
153  
154 <  void SelectionEvaluator::instructionDispatchLoop(OpenMDBitSet& bs, int frame){
154 >  void SelectionEvaluator::instructionDispatchLoop(SelectionSet& bs, int frame){
155      
156      while ( pc < aatoken.size()) {
157        statement = aatoken[pc++];
# Line 169 | Line 172 | namespace OpenMD {
172  
173    }
174  
175 <  OpenMDBitSet SelectionEvaluator::expression(const std::vector<Token>& code,
175 >  SelectionSet SelectionEvaluator::expression(const std::vector<Token>& code,
176                                               int pcStart) {
177 <    OpenMDBitSet bs;
178 <    std::stack<OpenMDBitSet> stack;
177 >    SelectionSet bs = createSelectionSets();
178 >    std::stack<SelectionSet> stack;
179 >    vector<int> bsSize = bs.size();
180    
181      for (unsigned int pc = pcStart; pc < code.size(); ++pc) {
182        Token instruction = code[pc];
# Line 184 | Line 188 | namespace OpenMD {
188          break;
189        case Token::all:
190          bs = allInstruction();
191 <        stack.push(bs);            
191 >        stack.push(bs);
192          break;
193        case Token::none:
194 <        bs = OpenMDBitSet(nStuntDouble);
195 <        stack.push(bs);            
194 >        bs = createSelectionSets();
195 >        stack.push(bs);
196          break;
197        case Token::opOr:
198 <        bs = stack.top();
199 <        stack.pop();
200 <        stack.top() |= bs;
198 >          bs= stack.top();
199 >          stack.pop();
200 >          stack.top() |= bs;
201          break;
202        case Token::opAnd:
203 <        bs = stack.top();
204 <        stack.pop();
205 <        stack.top() &= bs;
203 >          bs = stack.top();
204 >          stack.pop();
205 >          stack.top() &= bs;
206          break;
207        case Token::opNot:
208 <        stack.top().flip();
208 >          stack.top().flip();
209          break;
210        case Token::within:
211          withinInstruction(instruction, stack.top());
212          break;
213        case Token::hull:
214 <        stack.push(hull());
214 >          stack.push(hull());
215          break;
216          //case Token::selected:
217          //  stack.push(getSelectionSet());
218          //  break;
219        case Token::name:
220 <        stack.push(nameInstruction(boost::any_cast<std::string>(instruction.value)));
220 >          stack.push(nameInstruction(boost::any_cast<std::string>(instruction.value)));
221          break;
222        case Token::index:
223 <        stack.push(indexInstruction(instruction.value));
223 >          stack.push(indexInstruction(instruction.value));
224          break;
225        case Token::identifier:
226 <        stack.push(lookupValue(boost::any_cast<std::string>(instruction.value)));
226 >          stack.push(lookupValue(boost::any_cast<std::string>(instruction.value)));
227          break;
228        case Token::opLT:
229        case Token::opLE:
# Line 240 | Line 244 | namespace OpenMD {
244    }
245  
246  
247 <  OpenMDBitSet SelectionEvaluator::expression(const std::vector<Token>& code,
247 >  SelectionSet SelectionEvaluator::expression(const std::vector<Token>& code,
248                                                int pcStart, int frame) {
249 <    OpenMDBitSet bs;
250 <    std::stack<OpenMDBitSet> stack;
249 >    SelectionSet bs = createSelectionSets();
250 >    std::stack<SelectionSet> stack;
251    
252      for (unsigned int pc = pcStart; pc < code.size(); ++pc) {
253        Token instruction = code[pc];
# Line 258 | Line 262 | namespace OpenMD {
262          stack.push(bs);            
263          break;
264        case Token::none:
265 <        bs = OpenMDBitSet(nStuntDouble);
265 >        bs = SelectionSet(nObjects);
266          stack.push(bs);            
267          break;
268        case Token::opOr:
# Line 312 | Line 316 | namespace OpenMD {
316  
317  
318  
319 <  OpenMDBitSet SelectionEvaluator::comparatorInstruction(const Token& instruction) {
319 >  SelectionSet SelectionEvaluator::comparatorInstruction(const Token& instruction) {
320      int comparator = instruction.tok;
321      int property = instruction.intValue;
322      float comparisonValue = boost::any_cast<float>(instruction.value);
323 <    OpenMDBitSet bs(nStuntDouble);
323 >    SelectionSet bs = createSelectionSets();
324      bs.clearAll();
325      
326      SimInfo::MoleculeIterator mi;
# Line 329 | Line 333 | namespace OpenMD {
333      for (mol = info->beginMolecule(mi); mol != NULL;
334           mol = info->nextMolecule(mi)) {
335  
336 +      compareProperty(mol, bs, property, comparator, comparisonValue);
337 +
338        for(atom = mol->beginAtom(ai); atom != NULL; atom = mol->nextAtom(ai)) {
339          compareProperty(atom, bs, property, comparator, comparisonValue);
340        }
# Line 342 | Line 348 | namespace OpenMD {
348      return bs;
349    }
350  
351 <  OpenMDBitSet SelectionEvaluator::comparatorInstruction(const Token& instruction, int frame) {
351 >  SelectionSet SelectionEvaluator::comparatorInstruction(const Token& instruction, int frame) {
352      int comparator = instruction.tok;
353      int property = instruction.intValue;
354      float comparisonValue = boost::any_cast<float>(instruction.value);
355 <    OpenMDBitSet bs(nStuntDouble);
355 >    SelectionSet bs = createSelectionSets();
356      bs.clearAll();
357      
358      SimInfo::MoleculeIterator mi;
# Line 359 | Line 365 | namespace OpenMD {
365      for (mol = info->beginMolecule(mi); mol != NULL;
366           mol = info->nextMolecule(mi)) {
367  
368 +      compareProperty(mol, bs, property, comparator, comparisonValue, frame);
369 +
370        for(atom = mol->beginAtom(ai); atom != NULL; atom = mol->nextAtom(ai)) {
371          compareProperty(atom, bs, property, comparator, comparisonValue, frame);
372        }
# Line 372 | Line 380 | namespace OpenMD {
380      return bs;
381    }
382  
383 <  void SelectionEvaluator::compareProperty(StuntDouble* sd, OpenMDBitSet& bs,
383 >  void SelectionEvaluator::compareProperty(StuntDouble* sd, SelectionSet& bs,
384                                             int property, int comparator,
385                                             float comparisonValue) {
386      RealType propertyValue = 0.0;
387      Vector3d pos;
388 +
389      switch (property) {
390      case Token::mass:
391        propertyValue = sd->getMass();
# Line 446 | Line 455 | namespace OpenMD {
455        match = propertyValue != comparisonValue;
456        break;
457      }
458 +
459      if (match)
460 <      bs.setBitOn(sd->getGlobalIndex());
460 >      bs.bitsets_[STUNTDOUBLE].setBitOn(sd->getGlobalIndex());
461  
462      
463    }
464  
465 <  void SelectionEvaluator::compareProperty(StuntDouble* sd, OpenMDBitSet& bs,
465 >  void SelectionEvaluator::compareProperty(Molecule* mol, SelectionSet& bs,
466 >                                           int property, int comparator,
467 >                                           float comparisonValue) {
468 >    RealType propertyValue = 0.0;
469 >    Vector3d pos;
470 >
471 >    switch (property) {
472 >    case Token::mass:
473 >      propertyValue = mol->getMass();
474 >      break;
475 >    case Token::charge:
476 >      {
477 >        Molecule::AtomIterator ai;
478 >        Atom* atom;
479 >        for (atom = mol->beginAtom(ai); atom != NULL;
480 >             atom = mol->nextAtom(ai)) {
481 >          propertyValue+=  getCharge(atom);
482 >        }
483 >      }
484 >      break;
485 >    case Token::x:
486 >      propertyValue = mol->getCom().x();
487 >      break;
488 >    case Token::y:
489 >      propertyValue = mol->getCom().y();
490 >      break;
491 >    case Token::z:
492 >      propertyValue = mol->getCom().z();
493 >      break;
494 >    case Token::wrappedX:    
495 >      pos = mol->getCom();
496 >      info->getSnapshotManager()->getCurrentSnapshot()->wrapVector(pos);
497 >      propertyValue = pos.x();
498 >      break;
499 >    case Token::wrappedY:
500 >      pos = mol->getCom();
501 >      info->getSnapshotManager()->getCurrentSnapshot()->wrapVector(pos);
502 >      propertyValue = pos.y();
503 >      break;
504 >    case Token::wrappedZ:
505 >      pos = mol->getCom();
506 >      info->getSnapshotManager()->getCurrentSnapshot()->wrapVector(pos);
507 >      propertyValue = pos.z();
508 >      break;
509 >    case Token::r:
510 >      propertyValue = mol->getCom().length();
511 >      break;
512 >    default:
513 >      unrecognizedMoleculeProperty(property);
514 >    }
515 >    
516 >    bool match = false;
517 >    switch (comparator) {
518 >    case Token::opLT:
519 >      match = propertyValue < comparisonValue;
520 >      break;
521 >    case Token::opLE:
522 >      match = propertyValue <= comparisonValue;
523 >      break;
524 >    case Token::opGE:
525 >      match = propertyValue >= comparisonValue;
526 >      break;
527 >    case Token::opGT:
528 >      match = propertyValue > comparisonValue;
529 >      break;
530 >    case Token::opEQ:
531 >      match = propertyValue == comparisonValue;
532 >      break;
533 >    case Token::opNE:
534 >      match = propertyValue != comparisonValue;
535 >      break;
536 >    }
537 >    
538 >    if (match)
539 >      bs.bitsets_[MOLECULE].setBitOn(mol->getGlobalIndex());    
540 >  }
541 >
542 >  void SelectionEvaluator::compareProperty(Molecule* mol, SelectionSet& bs,
543                                             int property, int comparator,
544                                             float comparisonValue, int frame) {
545      RealType propertyValue = 0.0;
546      Vector3d pos;
547      switch (property) {
548      case Token::mass:
549 +      propertyValue = mol->getMass();
550 +      break;
551 +    case Token::charge:
552 +      {
553 +        Molecule::AtomIterator ai;
554 +        Atom* atom;
555 +        for (atom = mol->beginAtom(ai); atom != NULL;
556 +             atom = mol->nextAtom(ai)) {
557 +          propertyValue+=  getCharge(atom,frame);
558 +        }
559 +      }      
560 +      break;
561 +    case Token::x:
562 +      propertyValue = mol->getCom(frame).x();
563 +      break;
564 +    case Token::y:
565 +      propertyValue = mol->getCom(frame).y();
566 +      break;
567 +    case Token::z:
568 +      propertyValue = mol->getCom(frame).z();
569 +      break;
570 +    case Token::wrappedX:    
571 +      pos = mol->getCom(frame);
572 +      info->getSnapshotManager()->getSnapshot(frame)->wrapVector(pos);
573 +      propertyValue = pos.x();
574 +      break;
575 +    case Token::wrappedY:
576 +      pos = mol->getCom(frame);
577 +      info->getSnapshotManager()->getSnapshot(frame)->wrapVector(pos);
578 +      propertyValue = pos.y();
579 +      break;
580 +    case Token::wrappedZ:
581 +      pos = mol->getCom(frame);
582 +      info->getSnapshotManager()->getSnapshot(frame)->wrapVector(pos);
583 +      propertyValue = pos.z();
584 +      break;
585 +
586 +    case Token::r:
587 +      propertyValue = mol->getCom(frame).length();
588 +      break;
589 +    default:
590 +      unrecognizedMoleculeProperty(property);
591 +    }
592 +        
593 +    bool match = false;
594 +    switch (comparator) {
595 +    case Token::opLT:
596 +      match = propertyValue < comparisonValue;
597 +      break;
598 +    case Token::opLE:
599 +      match = propertyValue <= comparisonValue;
600 +      break;
601 +    case Token::opGE:
602 +      match = propertyValue >= comparisonValue;
603 +      break;
604 +    case Token::opGT:
605 +      match = propertyValue > comparisonValue;
606 +      break;
607 +    case Token::opEQ:
608 +      match = propertyValue == comparisonValue;
609 +      break;
610 +    case Token::opNE:
611 +      match = propertyValue != comparisonValue;
612 +      break;
613 +    }
614 +    if (match)
615 +      bs.bitsets_[MOLECULE].setBitOn(mol->getGlobalIndex());
616 +    
617 +    
618 +  }
619 +  void SelectionEvaluator::compareProperty(StuntDouble* sd, SelectionSet& bs,
620 +                                           int property, int comparator,
621 +                                           float comparisonValue, int frame) {
622 +    RealType propertyValue = 0.0;
623 +    Vector3d pos;
624 +    switch (property) {
625 +    case Token::mass:
626        propertyValue = sd->getMass();
627        break;
628      case Token::charge:
# Line 528 | Line 692 | namespace OpenMD {
692        break;
693      }
694      if (match)
695 <      bs.setBitOn(sd->getGlobalIndex());
532 <
533 <    
695 >      bs.bitsets_[STUNTDOUBLE].setBitOn(sd->getGlobalIndex());    
696    }
697 <
697 >  
698 >  
699    void SelectionEvaluator::withinInstruction(const Token& instruction,
700 <                                             OpenMDBitSet& bs){
700 >                                             SelectionSet& bs){
701      
702      boost::any withinSpec = instruction.value;
703      float distance;
# Line 551 | Line 714 | namespace OpenMD {
714    }
715  
716    void SelectionEvaluator::withinInstruction(const Token& instruction,
717 <                                             OpenMDBitSet& bs, int frame){
717 >                                             SelectionSet& bs, int frame){
718      
719      boost::any withinSpec = instruction.value;
720      float distance;
# Line 608 | Line 771 | namespace OpenMD {
771      }
772    }
773  
774 <  void SelectionEvaluator::select(OpenMDBitSet& bs){
774 >  void SelectionEvaluator::select(SelectionSet& bs){
775      bs = expression(statement, 1);
776    }
777  
778 <  void SelectionEvaluator::select(OpenMDBitSet& bs, int frame){
778 >  void SelectionEvaluator::select(SelectionSet& bs, int frame){
779      bs = expression(statement, 1, frame);
780    }
781    
782 <  OpenMDBitSet SelectionEvaluator::lookupValue(const std::string& variable){
782 >  SelectionSet SelectionEvaluator::lookupValue(const std::string& variable){
783      
784 <    OpenMDBitSet bs(nStuntDouble);
784 >    SelectionSet bs = createSelectionSets();
785      std::map<std::string, boost::any>::iterator i = variables.find(variable);
786      
787      if (i != variables.end()) {
788 <      if (i->second.type() == typeid(OpenMDBitSet)) {
789 <        return boost::any_cast<OpenMDBitSet>(i->second);
788 >      if (i->second.type() == typeid(SelectionSet)) {
789 >        return boost::any_cast<SelectionSet>(i->second);
790        } else if (i->second.type() ==  typeid(std::vector<Token>)){
791          bs = expression(boost::any_cast<std::vector<Token> >(i->second), 2);
792          i->second =  bs; /**@todo fixme */
# Line 636 | Line 799 | namespace OpenMD {
799      return bs;
800    }
801    
802 <  OpenMDBitSet SelectionEvaluator::nameInstruction(const std::string& name){    
802 >  SelectionSet SelectionEvaluator::nameInstruction(const std::string& name){    
803      return nameFinder.match(name);    
804    }    
805  
# Line 657 | Line 820 | namespace OpenMD {
820      //predefine();
821    }
822  
823 <  OpenMDBitSet SelectionEvaluator::evaluate() {
824 <    OpenMDBitSet bs(nStuntDouble);
823 >  SelectionSet SelectionEvaluator::createSelectionSets() {
824 >    SelectionSet ss(nObjects);
825 >    return ss;
826 >  }
827 >
828 >  SelectionSet SelectionEvaluator::evaluate() {
829 >    SelectionSet bs = createSelectionSets();
830      if (isLoaded_) {
831        pc = 0;
832        instructionDispatchLoop(bs);
# Line 666 | Line 834 | namespace OpenMD {
834      return bs;
835    }
836  
837 <  OpenMDBitSet SelectionEvaluator::evaluate(int frame) {
838 <    OpenMDBitSet bs(nStuntDouble);
837 >  SelectionSet SelectionEvaluator::evaluate(int frame) {
838 >    SelectionSet bs = createSelectionSets();
839      if (isLoaded_) {
840        pc = 0;
841        instructionDispatchLoop(bs, frame);
# Line 675 | Line 843 | namespace OpenMD {
843      return bs;
844    }
845  
846 <  OpenMDBitSet SelectionEvaluator::indexInstruction(const boost::any& value) {
847 <    OpenMDBitSet bs(nStuntDouble);
846 >  SelectionSet SelectionEvaluator::indexInstruction(const boost::any& value) {
847 >    SelectionSet bs = createSelectionSets();
848  
849      if (value.type() == typeid(int)) {
850        int index = boost::any_cast<int>(value);
851 <      if (index < 0 || index >= bs.size()) {
851 >      if (index < 0 || index >= bs.bitsets_[STUNTDOUBLE].size()) {
852          invalidIndex(index);
853        } else {
854          bs = indexFinder.find(index);
# Line 688 | Line 856 | namespace OpenMD {
856      } else if (value.type() == typeid(std::pair<int, int>)) {
857        std::pair<int, int> indexRange= boost::any_cast<std::pair<int, int> >(value);
858        assert(indexRange.first <= indexRange.second);
859 <      if (indexRange.first < 0 || indexRange.second >= bs.size()) {
859 >      if (indexRange.first < 0 ||
860 >          indexRange.second >= bs.bitsets_[STUNTDOUBLE].size()) {
861          invalidIndexRange(indexRange);
862        }else {
863          bs = indexFinder.find(indexRange.first, indexRange.second);
# Line 698 | Line 867 | namespace OpenMD {
867      return bs;
868    }
869  
870 <  OpenMDBitSet SelectionEvaluator::allInstruction() {
871 <    OpenMDBitSet bs(nStuntDouble);
870 >  SelectionSet SelectionEvaluator::allInstruction() {
871 >    SelectionSet ss = createSelectionSets();
872  
873      SimInfo::MoleculeIterator mi;
705 –    Molecule* mol;
874      Molecule::AtomIterator ai;
707 –    Atom* atom;
875      Molecule::RigidBodyIterator rbIter;
876 +    Molecule::BondIterator bondIter;
877 +    Molecule::BendIterator bendIter;
878 +    Molecule::TorsionIterator torsionIter;
879 +    Molecule::InversionIterator inversionIter;
880 +
881 +    Molecule* mol;
882 +    Atom* atom;
883      RigidBody* rb;
884 +    Bond* bond;
885 +    Bend* bend;
886 +    Torsion* torsion;
887 +    Inversion* inversion;    
888  
889      // Doing the loop insures that we're actually on this processor.
890  
891      for (mol = info->beginMolecule(mi); mol != NULL;
892           mol = info->nextMolecule(mi)) {
715 –
893        for(atom = mol->beginAtom(ai); atom != NULL; atom = mol->nextAtom(ai)) {
894 <        bs.setBitOn(atom->getGlobalIndex());
895 <      }
719 <    
894 >        ss.bitsets_[STUNTDOUBLE].setBitOn(atom->getGlobalIndex());
895 >      }    
896        for (rb = mol->beginRigidBody(rbIter); rb != NULL;
897             rb = mol->nextRigidBody(rbIter)) {
898 <        bs.setBitOn(rb->getGlobalIndex());
898 >        ss.bitsets_[STUNTDOUBLE].setBitOn(rb->getGlobalIndex());
899        }
900 +      for (bond = mol->beginBond(bondIter); bond != NULL;
901 +           bond = mol->nextBond(bondIter)) {
902 +        ss.bitsets_[BOND].setBitOn(bond->getGlobalIndex());
903 +      }  
904 +      for (bend = mol->beginBend(bendIter); bend != NULL;
905 +           bend = mol->nextBend(bendIter)) {
906 +        ss.bitsets_[BEND].setBitOn(bend->getGlobalIndex());
907 +      }  
908 +      for (torsion = mol->beginTorsion(torsionIter); torsion != NULL;
909 +           torsion = mol->nextTorsion(torsionIter)) {
910 +        ss.bitsets_[TORSION].setBitOn(torsion->getGlobalIndex());
911 +      }  
912 +      for (inversion = mol->beginInversion(inversionIter); inversion != NULL;
913 +           inversion = mol->nextInversion(inversionIter)) {
914 +        ss.bitsets_[INVERSION].setBitOn(inversion->getGlobalIndex());
915 +      }
916 +      ss.bitsets_[MOLECULE].setBitOn(mol->getGlobalIndex());
917      }
918  
919 <    return bs;
919 >    return ss;
920    }
921  
922 <  OpenMDBitSet SelectionEvaluator::hull() {
923 <    OpenMDBitSet bs(nStuntDouble);
922 >  SelectionSet SelectionEvaluator::hull() {
923 >    SelectionSet bs = createSelectionSets();
924      
925      bs = hullFinder.findHull();
926      surfaceArea_ = hullFinder.getSurfaceArea();
# Line 736 | Line 929 | namespace OpenMD {
929    }
930  
931  
932 <  OpenMDBitSet SelectionEvaluator::hull(int frame) {
933 <    OpenMDBitSet bs(nStuntDouble);
932 >  SelectionSet SelectionEvaluator::hull(int frame) {
933 >    SelectionSet bs = createSelectionSets();
934      
935      bs = hullFinder.findHull(frame);
936  

Diff Legend

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