Changes in src/analysis_bonds.cpp [68f03d:112b09]
- File:
-
- 1 edited
-
src/analysis_bonds.cpp (modified) (11 diffs)
Legend:
- Unmodified
- Added
- Removed
-
src/analysis_bonds.cpp
r68f03d r112b09 5 5 * Author: heber 6 6 */ 7 8 #include "Helpers/MemDebug.hpp" 7 9 8 10 #include "analysis_bonds.hpp" … … 26 28 Mean = 0.; 27 29 28 atom *Walker = mol->start;29 30 int AtomCount = 0; 30 while (Walker->next != mol->end) { 31 Walker = Walker->next; 32 const int count = Walker->ListOfBonds.size(); 31 for (molecule::const_iterator iter = mol->begin(); iter != mol->end(); ++iter) { 32 const int count = (*iter)->ListOfBonds.size(); 33 33 if (Max < count) 34 34 Max = count; … … 51 51 * \param &Max maximum distance on return, 0 if no bond between the two elements 52 52 */ 53 void MinMeanMaxBondDistanceBetweenElements(const molecule *mol, element *type1,element *type2, double &Min, double &Mean, double &Max)53 void MinMeanMaxBondDistanceBetweenElements(const molecule *mol, const element *type1, const element *type2, double &Min, double &Mean, double &Max) 54 54 { 55 55 Min = 2e+6; … … 58 58 59 59 int AtomNo = 0; 60 atom *Walker = mol->start; 61 while (Walker->next != mol->end) { 62 Walker = Walker->next; 63 if (Walker->type == type1) 64 for (BondList::const_iterator BondRunner = Walker->ListOfBonds.begin(); BondRunner != Walker->ListOfBonds.end(); BondRunner++) 65 if ((*BondRunner)->GetOtherAtom(Walker)->type == type2) { 60 for (molecule::const_iterator iter = mol->begin(); iter != mol->end(); ++iter) { 61 if ((*iter)->type == type1) 62 for (BondList::const_iterator BondRunner = (*iter)->ListOfBonds.begin(); BondRunner != (*iter)->ListOfBonds.end(); BondRunner++) 63 if ((*BondRunner)->GetOtherAtom((*iter))->type == type2) { 66 64 const double distance = (*BondRunner)->GetDistanceSquared(); 67 65 if (Min > distance) … … 124 122 * \param *InterfaceElement or NULL 125 123 */ 126 int CountHydrogenBridgeBonds(MoleculeListClass *molecules, element * InterfaceElement = NULL) 127 { 128 atom *Walker = NULL; 129 atom *Runner = NULL; 124 int CountHydrogenBridgeBonds(MoleculeListClass *molecules, const element * InterfaceElement = NULL) 125 { 130 126 int count = 0; 131 127 int OtherHydrogens = 0; … … 133 129 bool InterfaceFlag = false; 134 130 bool OtherHydrogenFlag = true; 135 for (MoleculeList::const_iterator MolWalker = molecules->ListOfMolecules.begin();MolWalker != molecules->ListOfMolecules.end(); MolWalker++) { 136 Walker = (*MolWalker)->start; 137 while (Walker->next != (*MolWalker)->end) { 138 Walker = Walker->next; 139 for (MoleculeList::const_iterator MolRunner = molecules->ListOfMolecules.begin();MolRunner != molecules->ListOfMolecules.end(); MolRunner++) { 140 Runner = (*MolRunner)->start; 141 while (Runner->next != (*MolRunner)->end) { 142 Runner = Runner->next; 143 if ((Walker->type->Z == 8) && (Runner->type->Z == 8)) { 131 for (MoleculeList::const_iterator MolWalker = molecules->ListOfMolecules.begin();MolWalker != molecules->ListOfMolecules.end(); ++MolWalker) { 132 molecule::iterator Walker = (*MolWalker)->begin(); 133 for(;Walker!=(*MolWalker)->end();++Walker){ 134 for (MoleculeList::const_iterator MolRunner = molecules->ListOfMolecules.begin();MolRunner != molecules->ListOfMolecules.end(); ++MolRunner) { 135 molecule::iterator Runner = (*MolRunner)->begin(); 136 for(;Runner!=(*MolRunner)->end();++Runner){ 137 if (((*Walker)->type->Z == 8) && ((*Runner)->type->Z == 8)) { 144 138 // check distance 145 const double distance = Runner->x.DistanceSquared(Walker->x);139 const double distance = (*Runner)->x.DistanceSquared((*Walker)->x); 146 140 if ((distance > MYEPSILON) && (distance < HBRIDGEDISTANCE*HBRIDGEDISTANCE)) { // distance >0 means different atoms 147 141 // on other atom(Runner) we check for bond to interface element and … … 151 145 OtherHydrogens = 0; 152 146 InterfaceFlag = (InterfaceElement == NULL); 153 for (BondList::const_iterator BondRunner = Runner->ListOfBonds.begin(); BondRunner != Runner->ListOfBonds.end(); BondRunner++) {154 atom * const OtherAtom = (*BondRunner)->GetOtherAtom( Runner);147 for (BondList::const_iterator BondRunner = (*Runner)->ListOfBonds.begin(); BondRunner != (*Runner)->ListOfBonds.end(); BondRunner++) { 148 atom * const OtherAtom = (*BondRunner)->GetOtherAtom(*Runner); 155 149 // if hydrogen, check angle to be greater(!) than 30 degrees 156 150 if (OtherAtom->type->Z == 1) { 157 const double angle = CalculateAngle(&OtherAtom->x, & Runner->x, &Walker->x);151 const double angle = CalculateAngle(&OtherAtom->x, &(*Runner)->x, &(*Walker)->x); 158 152 OtherHydrogenFlag = OtherHydrogenFlag && (angle > M_PI*(30./180.) + MYEPSILON); 159 153 Otherangle += angle; … … 176 170 if (InterfaceFlag && OtherHydrogenFlag) { 177 171 // on this element (Walker) we check for bond to hydrogen, i.e. part of water molecule 178 for (BondList::const_iterator BondRunner = Walker->ListOfBonds.begin(); BondRunner != Walker->ListOfBonds.end(); BondRunner++) {179 atom * const OtherAtom = (*BondRunner)->GetOtherAtom( Walker);172 for (BondList::const_iterator BondRunner = (*Walker)->ListOfBonds.begin(); BondRunner != (*Walker)->ListOfBonds.end(); BondRunner++) { 173 atom * const OtherAtom = (*BondRunner)->GetOtherAtom(*Walker); 180 174 if (OtherAtom->type->Z == 1) { 181 175 // check angle 182 if (CheckHydrogenBridgeBondAngle( Walker, OtherAtom,Runner)) {183 DoLog(1) && (Log() << Verbose(1) << Walker->getName() << ", " << OtherAtom->getName() << " and " << Runner->getName() << " has a hydrogen bridge bond with distance " << sqrt(distance) << " and angle " << CalculateAngle(&OtherAtom->x, &Walker->x, &Runner->x)*(180./M_PI) << "." << endl);176 if (CheckHydrogenBridgeBondAngle(*Walker, OtherAtom, *Runner)) { 177 DoLog(1) && (Log() << Verbose(1) << (*Walker)->getName() << ", " << OtherAtom->getName() << " and " << (*Runner)->getName() << " has a hydrogen bridge bond with distance " << sqrt(distance) << " and angle " << CalculateAngle(&OtherAtom->x, &(*Walker)->x, &(*Runner)->x)*(180./M_PI) << "." << endl); 184 178 count++; 185 179 break; … … 205 199 int CountBondsOfTwo(MoleculeListClass * const molecules, const element * const first, const element * const second) 206 200 { 207 atom *Walker = NULL;208 201 int count = 0; 209 202 210 203 for (MoleculeList::const_iterator MolWalker = molecules->ListOfMolecules.begin();MolWalker != molecules->ListOfMolecules.end(); MolWalker++) { 211 Walker = (*MolWalker)->start;212 while (Walker->next != (*MolWalker)->end){213 Walker = Walker->next;214 if (( Walker->type == first) || (Walker->type == second)) { // first element matches215 for (BondList::const_iterator BondRunner = Walker->ListOfBonds.begin(); BondRunner != Walker->ListOfBonds.end(); BondRunner++) {216 atom * const OtherAtom = (*BondRunner)->GetOtherAtom( Walker);217 if (((OtherAtom->type == first) || (OtherAtom->type == second)) && ( Walker->nr < OtherAtom->nr)) {204 molecule::iterator Walker = (*MolWalker)->begin(); 205 for(;Walker!=(*MolWalker)->end();++Walker){ 206 atom * theAtom = *Walker; 207 if ((theAtom->type == first) || (theAtom->type == second)) { // first element matches 208 for (BondList::const_iterator BondRunner = theAtom->ListOfBonds.begin(); BondRunner != theAtom->ListOfBonds.end(); BondRunner++) { 209 atom * const OtherAtom = (*BondRunner)->GetOtherAtom(theAtom); 210 if (((OtherAtom->type == first) || (OtherAtom->type == second)) && (theAtom->nr < OtherAtom->nr)) { 218 211 count++; 219 212 DoLog(1) && (Log() << Verbose(1) << first->name << "-" << second->name << " bond found between " << *Walker << " and " << *OtherAtom << "." << endl); … … 240 233 bool MatchFlag[2]; 241 234 bool result = false; 242 atom *Walker = NULL;243 235 const element * ElementArray[2]; 244 236 ElementArray[0] = first; … … 246 238 247 239 for (MoleculeList::const_iterator MolWalker = molecules->ListOfMolecules.begin();MolWalker != molecules->ListOfMolecules.end(); MolWalker++) { 248 Walker = (*MolWalker)->start;249 while (Walker->next != (*MolWalker)->end){250 Walker = Walker->next;251 if ( Walker->type == second) { // first element matches240 molecule::iterator Walker = (*MolWalker)->begin(); 241 for(;Walker!=(*MolWalker)->end();++Walker){ 242 atom *theAtom = *Walker; 243 if (theAtom->type == second) { // first element matches 252 244 for (int i=0;i<2;i++) 253 245 MatchFlag[i] = false; 254 for (BondList::const_iterator BondRunner = Walker->ListOfBonds.begin(); BondRunner != Walker->ListOfBonds.end(); BondRunner++) {255 atom * const OtherAtom = (*BondRunner)->GetOtherAtom( Walker);246 for (BondList::const_iterator BondRunner = theAtom->ListOfBonds.begin(); BondRunner != theAtom->ListOfBonds.end(); BondRunner++) { 247 atom * const OtherAtom = (*BondRunner)->GetOtherAtom(theAtom); 256 248 for (int i=0;i<2;i++) 257 249 if ((!MatchFlag[i]) && (OtherAtom->type == ElementArray[i])) {
Note:
See TracChangeset
for help on using the changeset viewer.
