| 1 | /* | 
|---|
| 2 | * Project: MoleCuilder | 
|---|
| 3 | * Description: creates and alters molecular systems | 
|---|
| 4 | * Copyright (C)  2010-2012 University of Bonn. All rights reserved. | 
|---|
| 5 | * Please see the LICENSE file or "Copyright notice" in builder.cpp for details. | 
|---|
| 6 | */ | 
|---|
| 7 |  | 
|---|
| 8 | /* | 
|---|
| 9 | * TremoloParser.cpp | 
|---|
| 10 | * | 
|---|
| 11 | *  Created on: Mar 2, 2010 | 
|---|
| 12 | *      Author: metzler | 
|---|
| 13 | */ | 
|---|
| 14 |  | 
|---|
| 15 | // include config.h | 
|---|
| 16 | #ifdef HAVE_CONFIG_H | 
|---|
| 17 | #include <config.h> | 
|---|
| 18 | #endif | 
|---|
| 19 |  | 
|---|
| 20 | #include "CodePatterns/MemDebug.hpp" | 
|---|
| 21 |  | 
|---|
| 22 | #include "CodePatterns/Assert.hpp" | 
|---|
| 23 | #include "CodePatterns/Log.hpp" | 
|---|
| 24 | #include "CodePatterns/toString.hpp" | 
|---|
| 25 | #include "CodePatterns/Verbose.hpp" | 
|---|
| 26 |  | 
|---|
| 27 | #include "TremoloParser.hpp" | 
|---|
| 28 |  | 
|---|
| 29 | #include "Atom/atom.hpp" | 
|---|
| 30 | #include "Bond/bond.hpp" | 
|---|
| 31 | #include "Box.hpp" | 
|---|
| 32 | #include "Descriptors/AtomIdDescriptor.hpp" | 
|---|
| 33 | #include "Element/element.hpp" | 
|---|
| 34 | #include "Element/periodentafel.hpp" | 
|---|
| 35 | #include "LinearAlgebra/RealSpaceMatrix.hpp" | 
|---|
| 36 | #include "molecule.hpp" | 
|---|
| 37 | #include "MoleculeListClass.hpp" | 
|---|
| 38 | #include "World.hpp" | 
|---|
| 39 | #include "WorldTime.hpp" | 
|---|
| 40 |  | 
|---|
| 41 |  | 
|---|
| 42 | #include <algorithm> | 
|---|
| 43 | #include <boost/lexical_cast.hpp> | 
|---|
| 44 | #include <boost/tokenizer.hpp> | 
|---|
| 45 | #include <iostream> | 
|---|
| 46 | #include <iomanip> | 
|---|
| 47 | #include <map> | 
|---|
| 48 | #include <sstream> | 
|---|
| 49 | #include <vector> | 
|---|
| 50 |  | 
|---|
| 51 | // declare specialized static variables | 
|---|
| 52 | const std::string FormatParserTrait<tremolo>::name = "tremolo"; | 
|---|
| 53 | const std::string FormatParserTrait<tremolo>::suffix = "data"; | 
|---|
| 54 | const ParserTypes FormatParserTrait<tremolo>::type = tremolo; | 
|---|
| 55 |  | 
|---|
| 56 | /** | 
|---|
| 57 | * Constructor. | 
|---|
| 58 | */ | 
|---|
| 59 | FormatParser< tremolo >::FormatParser()  : | 
|---|
| 60 | FormatParser_common(NULL) | 
|---|
| 61 | { | 
|---|
| 62 | knownKeys["x"] = TremoloKey::x; | 
|---|
| 63 | knownKeys["u"] = TremoloKey::u; | 
|---|
| 64 | knownKeys["F"] = TremoloKey::F; | 
|---|
| 65 | knownKeys["stress"] = TremoloKey::stress; | 
|---|
| 66 | knownKeys["Id"] = TremoloKey::Id; | 
|---|
| 67 | knownKeys["neighbors"] = TremoloKey::neighbors; | 
|---|
| 68 | knownKeys["imprData"] = TremoloKey::imprData; | 
|---|
| 69 | knownKeys["GroupMeasureTypeNo"] = TremoloKey::GroupMeasureTypeNo; | 
|---|
| 70 | knownKeys["type"] = TremoloKey::type; | 
|---|
| 71 | knownKeys["extType"] = TremoloKey::extType; | 
|---|
| 72 | knownKeys["name"] = TremoloKey::name; | 
|---|
| 73 | knownKeys["resName"] = TremoloKey::resName; | 
|---|
| 74 | knownKeys["chainID"] = TremoloKey::chainID; | 
|---|
| 75 | knownKeys["resSeq"] = TremoloKey::resSeq; | 
|---|
| 76 | knownKeys["occupancy"] = TremoloKey::occupancy; | 
|---|
| 77 | knownKeys["tempFactor"] = TremoloKey::tempFactor; | 
|---|
| 78 | knownKeys["segID"] = TremoloKey::segID; | 
|---|
| 79 | knownKeys["Charge"] = TremoloKey::Charge; | 
|---|
| 80 | knownKeys["charge"] = TremoloKey::charge; | 
|---|
| 81 | knownKeys["GrpTypeNo"] = TremoloKey::GrpTypeNo; | 
|---|
| 82 | knownKeys["torsion"] = TremoloKey::torsion; | 
|---|
| 83 | knownKeys[" "] = TremoloKey::noKey; // with this we can detect invalid keys | 
|---|
| 84 |  | 
|---|
| 85 | createKnownTypesByIdentity(); | 
|---|
| 86 |  | 
|---|
| 87 | // invert knownKeys for debug output | 
|---|
| 88 | for (std::map<std::string, TremoloKey::atomDataKey>::iterator iter = knownKeys.begin(); iter != knownKeys.end(); ++iter) | 
|---|
| 89 | knownKeyNames.insert( make_pair( iter->second, iter->first) ); | 
|---|
| 90 |  | 
|---|
| 91 | additionalAtomData.clear(); | 
|---|
| 92 | } | 
|---|
| 93 |  | 
|---|
| 94 | /** | 
|---|
| 95 | * Destructor. | 
|---|
| 96 | */ | 
|---|
| 97 | FormatParser< tremolo >::~FormatParser() | 
|---|
| 98 | { | 
|---|
| 99 | usedFields_save.clear(); | 
|---|
| 100 | additionalAtomData.clear(); | 
|---|
| 101 | knownKeys.clear(); | 
|---|
| 102 | } | 
|---|
| 103 |  | 
|---|
| 104 | /** | 
|---|
| 105 | * Loads atoms from a tremolo-formatted file. | 
|---|
| 106 | * | 
|---|
| 107 | * \param tremolo file | 
|---|
| 108 | */ | 
|---|
| 109 | void FormatParser< tremolo >::load(istream* file) { | 
|---|
| 110 | std::string line; | 
|---|
| 111 | std::string::size_type location; | 
|---|
| 112 |  | 
|---|
| 113 | // reset the id maps | 
|---|
| 114 | resetIdAssociations(); | 
|---|
| 115 |  | 
|---|
| 116 | molecule *newmol = World::getInstance().createMolecule(); | 
|---|
| 117 | newmol->ActiveFlag = true; | 
|---|
| 118 | // TODO: Remove the insertion into molecule when saving does not depend on them anymore. Also, remove molecule.hpp include | 
|---|
| 119 | World::getInstance().getMolecules()->insert(newmol); | 
|---|
| 120 | while (file->good()) { | 
|---|
| 121 | std::getline(*file, line, '\n'); | 
|---|
| 122 | // we only parse in the first ATOMDATA line | 
|---|
| 123 | if (usedFields_load.empty()) { | 
|---|
| 124 | location = line.find("ATOMDATA", 0); | 
|---|
| 125 | if (location != string::npos) { | 
|---|
| 126 | parseAtomDataKeysLine(line, location + 8, usedFields_load); | 
|---|
| 127 | } | 
|---|
| 128 | } | 
|---|
| 129 | if (line.length() > 0 && line.at(0) != '#') { | 
|---|
| 130 | readAtomDataLine(line, newmol); | 
|---|
| 131 | } | 
|---|
| 132 | } | 
|---|
| 133 | LOG(3, "DEBUG: Local usedFields is: " << usedFields_load); | 
|---|
| 134 |  | 
|---|
| 135 | // refresh atom::nr and atom::name | 
|---|
| 136 | std::vector<atomId_t> atoms(newmol->getAtomCount()); | 
|---|
| 137 | std::transform(newmol->begin(), newmol->end(), atoms.begin(), | 
|---|
| 138 | boost::bind(&atom::getId, _1)); | 
|---|
| 139 | processNeighborInformation(atoms); | 
|---|
| 140 | adaptImprData(); | 
|---|
| 141 | adaptTorsion(); | 
|---|
| 142 |  | 
|---|
| 143 | // append usedFields to global usedFields, is made unique on save, clear after use | 
|---|
| 144 | usedFields_save.insert(usedFields_save.end(), usedFields_load.begin(), usedFields_load.end()); | 
|---|
| 145 | usedFields_load.clear(); | 
|---|
| 146 | } | 
|---|
| 147 |  | 
|---|
| 148 | /** | 
|---|
| 149 | * Saves the \a atoms into as a tremolo file. | 
|---|
| 150 | * | 
|---|
| 151 | * \param file where to save the state | 
|---|
| 152 | * \param atoms atoms to store | 
|---|
| 153 | */ | 
|---|
| 154 | void FormatParser< tremolo >::save(std::ostream* file, const std::vector<atom *> &AtomList) { | 
|---|
| 155 | LOG(0, "Saving changes to tremolo."); | 
|---|
| 156 |  | 
|---|
| 157 | // install default usedFields if empty so far | 
|---|
| 158 | if (usedFields_save.empty()) { | 
|---|
| 159 | // default behavior: use all possible keys on output | 
|---|
| 160 | for (std::map<std::string, TremoloKey::atomDataKey>::iterator iter = knownKeys.begin(); | 
|---|
| 161 | iter != knownKeys.end(); ++iter) | 
|---|
| 162 | if (iter->second != TremoloKey::noKey) // don't add noKey | 
|---|
| 163 | usedFields_save.push_back(iter->first); | 
|---|
| 164 | } | 
|---|
| 165 | // make present usedFields_save unique | 
|---|
| 166 | makeUsedFieldsUnique(usedFields_save); | 
|---|
| 167 | LOG(1, "DEBUG: Global (with unique entries) usedFields_save is: " << usedFields_save); | 
|---|
| 168 |  | 
|---|
| 169 | // distribute ids continuously | 
|---|
| 170 | distributeContinuousIds(AtomList); | 
|---|
| 171 |  | 
|---|
| 172 | // store atomdata | 
|---|
| 173 | save_AtomDataLine(file); | 
|---|
| 174 |  | 
|---|
| 175 | // store box | 
|---|
| 176 | save_BoxLine(file); | 
|---|
| 177 |  | 
|---|
| 178 | // store particles | 
|---|
| 179 | for (std::vector<atom*>::const_iterator  atomIt = AtomList.begin(); | 
|---|
| 180 | atomIt != AtomList.end(); ++atomIt) | 
|---|
| 181 | saveLine(file, *atomIt); | 
|---|
| 182 | } | 
|---|
| 183 |  | 
|---|
| 184 | /** Helper function to make \given fields unique while preserving the order of first appearance. | 
|---|
| 185 | * | 
|---|
| 186 | * As std::unique only removes element if equal to predecessor, a vector is only | 
|---|
| 187 | * made unique if sorted beforehand. But sorting would destroy order of first | 
|---|
| 188 | * appearance, hence we do the sorting on a temporary field and add the unique | 
|---|
| 189 | * elements in the order as in \a fields. | 
|---|
| 190 | * | 
|---|
| 191 | * @param fields usedFields to make unique while preserving order of appearance | 
|---|
| 192 | */ | 
|---|
| 193 | void FormatParser< tremolo >::makeUsedFieldsUnique(usedFields_t &fields) | 
|---|
| 194 | { | 
|---|
| 195 | // std::unique only removes if predecessor is equal, not over whole range, hence do it manually | 
|---|
| 196 | usedFields_t temp_fields(usedFields_save); | 
|---|
| 197 | std::sort(temp_fields.begin(), temp_fields.end()); | 
|---|
| 198 | usedFields_t::iterator it = | 
|---|
| 199 | std::unique(temp_fields.begin(), temp_fields.end()); // skips all duplicates in the vector | 
|---|
| 200 | temp_fields.erase(it, temp_fields.end()); | 
|---|
| 201 | usedFields_t usedfields(usedFields_save); | 
|---|
| 202 | usedFields_save.clear(); | 
|---|
| 203 | usedFields_save.reserve(temp_fields.size()); | 
|---|
| 204 | // now go through each usedFields entry, check if in temp_fields and remove there on first occurence | 
|---|
| 205 | for (usedFields_t::const_iterator iter = usedfields.begin(); | 
|---|
| 206 | iter != usedfields.end(); ++iter) { | 
|---|
| 207 | usedFields_t::iterator uniqueiter = | 
|---|
| 208 | std::find(temp_fields.begin(), temp_fields.end(), *iter); | 
|---|
| 209 | if (uniqueiter != temp_fields.end()) { | 
|---|
| 210 | usedFields_save.push_back(*iter); | 
|---|
| 211 | // add only once to ATOMDATA | 
|---|
| 212 | temp_fields.erase(uniqueiter); | 
|---|
| 213 | } | 
|---|
| 214 | } | 
|---|
| 215 | ASSERT( temp_fields.empty(), | 
|---|
| 216 | "FormatParser< tremolo >::save() - still unique entries left in temp_fields after unique?"); | 
|---|
| 217 | } | 
|---|
| 218 |  | 
|---|
| 219 |  | 
|---|
| 220 | /** Resets and distributes the indices continuously. | 
|---|
| 221 | * | 
|---|
| 222 | * \param atoms atoms to store | 
|---|
| 223 | */ | 
|---|
| 224 | void FormatParser< tremolo >::distributeContinuousIds(const std::vector<atom *> &AtomList) | 
|---|
| 225 | { | 
|---|
| 226 | resetIdAssociations(); | 
|---|
| 227 | atomId_t lastid = 0; | 
|---|
| 228 | for (std::vector<atom*>::const_iterator atomIt = AtomList.begin(); | 
|---|
| 229 | atomIt != AtomList.end(); ++atomIt) | 
|---|
| 230 | associateLocaltoGlobalId(++lastid, (*atomIt)->getId()); | 
|---|
| 231 | } | 
|---|
| 232 |  | 
|---|
| 233 | /** Store Atomdata line to \a file. | 
|---|
| 234 | * | 
|---|
| 235 | * @param file output stream | 
|---|
| 236 | */ | 
|---|
| 237 | void FormatParser< tremolo >::save_AtomDataLine(std::ostream* file) const | 
|---|
| 238 | { | 
|---|
| 239 | *file << "# ATOMDATA"; | 
|---|
| 240 | for (usedFields_t::const_iterator it=usedFields_save.begin(); | 
|---|
| 241 | it != usedFields_save.end(); ++it) | 
|---|
| 242 | *file << "\t" << *it; | 
|---|
| 243 | *file << std::endl; | 
|---|
| 244 | } | 
|---|
| 245 |  | 
|---|
| 246 | /** Store Box info to \a file | 
|---|
| 247 | * | 
|---|
| 248 | * @param file output stream | 
|---|
| 249 | */ | 
|---|
| 250 | void FormatParser< tremolo >::save_BoxLine(std::ostream* file) const | 
|---|
| 251 | { | 
|---|
| 252 | *file << "# Box"; | 
|---|
| 253 | const RealSpaceMatrix &M = World::getInstance().getDomain().getM(); | 
|---|
| 254 | for (size_t i=0; i<NDIM;++i) | 
|---|
| 255 | for (size_t j=0; j<NDIM;++j) | 
|---|
| 256 | *file << "\t" << M.at(i,j); | 
|---|
| 257 | *file << std::endl; | 
|---|
| 258 | } | 
|---|
| 259 |  | 
|---|
| 260 | /** Add default info, when new atom is added to World. | 
|---|
| 261 | * | 
|---|
| 262 | * @param id of atom | 
|---|
| 263 | */ | 
|---|
| 264 | void FormatParser< tremolo >::AtomInserted(atomId_t id) | 
|---|
| 265 | { | 
|---|
| 266 | std::map<const atomId_t, TremoloAtomInfoContainer>::iterator iter = additionalAtomData.find(id); | 
|---|
| 267 | ASSERT(iter == additionalAtomData.end(), | 
|---|
| 268 | "FormatParser< tremolo >::AtomInserted() - additionalAtomData already present for newly added atom " | 
|---|
| 269 | +toString(id)+"."); | 
|---|
| 270 | // don't add entry, as this gives a default resSeq of 0 not the molecule id | 
|---|
| 271 | // additionalAtomData.insert( std::make_pair(id, TremoloAtomInfoContainer()) ); | 
|---|
| 272 | } | 
|---|
| 273 |  | 
|---|
| 274 | /** Remove additional AtomData info, when atom has been removed from World. | 
|---|
| 275 | * | 
|---|
| 276 | * @param id of atom | 
|---|
| 277 | */ | 
|---|
| 278 | void FormatParser< tremolo >::AtomRemoved(atomId_t id) | 
|---|
| 279 | { | 
|---|
| 280 | std::map<const atomId_t, TremoloAtomInfoContainer>::iterator iter = additionalAtomData.find(id); | 
|---|
| 281 | // as we do not insert AtomData on AtomInserted, we cannot be assured of its presence | 
|---|
| 282 | //  ASSERT(iter != additionalAtomData.end(), | 
|---|
| 283 | //      "FormatParser< tremolo >::AtomRemoved() - additionalAtomData is not present for atom " | 
|---|
| 284 | //      +toString(id)+" to remove."); | 
|---|
| 285 | if (iter != additionalAtomData.end()) | 
|---|
| 286 | additionalAtomData.erase(iter); | 
|---|
| 287 | } | 
|---|
| 288 |  | 
|---|
| 289 | /** | 
|---|
| 290 | * Writes one line of tremolo-formatted data to the provided stream. | 
|---|
| 291 | * | 
|---|
| 292 | * \param stream where to write the line to | 
|---|
| 293 | * \param reference to the atom of which information should be written | 
|---|
| 294 | */ | 
|---|
| 295 | void FormatParser< tremolo >::saveLine(std::ostream* file, const atom* currentAtom) | 
|---|
| 296 | { | 
|---|
| 297 | TremoloKey::atomDataKey currentField; | 
|---|
| 298 |  | 
|---|
| 299 | LOG(4, "INFO: Saving atom " << *currentAtom << ", its father id is " << currentAtom->GetTrueFather()->getId()); | 
|---|
| 300 |  | 
|---|
| 301 | for (usedFields_t::iterator it = usedFields_save.begin(); it != usedFields_save.end(); it++) { | 
|---|
| 302 | currentField = knownKeys[it->substr(0, it->find("="))]; | 
|---|
| 303 | switch (currentField) { | 
|---|
| 304 | case TremoloKey::x : | 
|---|
| 305 | // for the moment, assume there are always three dimensions | 
|---|
| 306 | LOG(3, "Writing for type " << knownKeyNames[currentField] << ": " << currentAtom->getPosition()); | 
|---|
| 307 | *file << currentAtom->at(0) << "\t"; | 
|---|
| 308 | *file << currentAtom->at(1) << "\t"; | 
|---|
| 309 | *file << currentAtom->at(2) << "\t"; | 
|---|
| 310 | break; | 
|---|
| 311 | case TremoloKey::u : | 
|---|
| 312 | // for the moment, assume there are always three dimensions | 
|---|
| 313 | LOG(3, "Writing for type " << knownKeyNames[currentField] << ": " << currentAtom->getAtomicVelocity()); | 
|---|
| 314 | *file << currentAtom->getAtomicVelocity()[0] << "\t"; | 
|---|
| 315 | *file << currentAtom->getAtomicVelocity()[1] << "\t"; | 
|---|
| 316 | *file << currentAtom->getAtomicVelocity()[2] << "\t"; | 
|---|
| 317 | break; | 
|---|
| 318 | case TremoloKey::type : | 
|---|
| 319 | if (additionalAtomData.count(currentAtom->getId())) { | 
|---|
| 320 | if (additionalAtomData[currentAtom->getId()].get(currentField) != "-") { | 
|---|
| 321 | LOG(3, "Writing for type " << knownKeyNames[currentField] << ": " << additionalAtomData[currentAtom->getId()].get(currentField)); | 
|---|
| 322 | *file << additionalAtomData[currentAtom->getId()].get(currentField) << "\t"; | 
|---|
| 323 | } else { | 
|---|
| 324 | LOG(3, "Writing for type " << knownKeyNames[currentField] << " default value: " << currentAtom->getType()->getSymbol()); | 
|---|
| 325 | *file << currentAtom->getType()->getSymbol() << "\t"; | 
|---|
| 326 | } | 
|---|
| 327 | } else if (additionalAtomData.count(currentAtom->GetTrueFather()->getId())) { | 
|---|
| 328 | if (additionalAtomData[currentAtom->GetTrueFather()->getId()].get(currentField) != "-") { | 
|---|
| 329 | LOG(3, "Writing for type " << knownKeyNames[currentField] << " stuff from father: " << additionalAtomData[currentAtom->GetTrueFather()->getId()].get(currentField)); | 
|---|
| 330 | *file << additionalAtomData[currentAtom->GetTrueFather()->getId()].get(currentField) << "\t"; | 
|---|
| 331 | } else { | 
|---|
| 332 | LOG(3, "Writing for type " << knownKeyNames[currentField] << " default value from father: " << currentAtom->GetTrueFather()->getType()->getSymbol()); | 
|---|
| 333 | *file << currentAtom->GetTrueFather()->getType()->getSymbol() << "\t"; | 
|---|
| 334 | } | 
|---|
| 335 | } else { | 
|---|
| 336 | LOG(3, "Writing for type " << knownKeyNames[currentField] << " its default value: " << currentAtom->getType()->getSymbol()); | 
|---|
| 337 | *file << currentAtom->getType()->getSymbol() << "\t"; | 
|---|
| 338 | } | 
|---|
| 339 | break; | 
|---|
| 340 | case TremoloKey::Id : | 
|---|
| 341 | LOG(3, "Writing for type " << knownKeyNames[currentField] << ": " << currentAtom->getId()+1); | 
|---|
| 342 | *file << getLocalId(currentAtom->getId()) << "\t"; | 
|---|
| 343 | break; | 
|---|
| 344 | case TremoloKey::neighbors : | 
|---|
| 345 | LOG(3, "Writing type " << knownKeyNames[currentField]); | 
|---|
| 346 | writeNeighbors(file, atoi(it->substr(it->find("=") + 1, 1).c_str()), currentAtom); | 
|---|
| 347 | break; | 
|---|
| 348 | case TremoloKey::resSeq : | 
|---|
| 349 | if (additionalAtomData.count(currentAtom->getId())) { | 
|---|
| 350 | LOG(3, "Writing for type " << knownKeyNames[currentField] << ": " << additionalAtomData[currentAtom->getId()].get(currentField)); | 
|---|
| 351 | *file << additionalAtomData[currentAtom->getId()].get(currentField); | 
|---|
| 352 | } else if (currentAtom->getMolecule() != NULL) { | 
|---|
| 353 | LOG(3, "Writing for type " << knownKeyNames[currentField] << " its own id: " << currentAtom->getMolecule()->getId()+1); | 
|---|
| 354 | *file << setw(4) << currentAtom->getMolecule()->getId()+1; | 
|---|
| 355 | } else { | 
|---|
| 356 | LOG(3, "Writing for type " << knownKeyNames[currentField] << " default value: " << defaultAdditionalData.get(currentField)); | 
|---|
| 357 | *file << defaultAdditionalData.get(currentField); | 
|---|
| 358 | } | 
|---|
| 359 | *file << "\t"; | 
|---|
| 360 | break; | 
|---|
| 361 | case TremoloKey::charge : | 
|---|
| 362 | if (currentAtom->getCharge() == 0.) { | 
|---|
| 363 | if (additionalAtomData.count(currentAtom->getId())) { | 
|---|
| 364 | LOG(3, "Writing for type " << knownKeyNames[currentField] << ": " << additionalAtomData[currentAtom->getId()].get(currentField)); | 
|---|
| 365 | *file << additionalAtomData[currentAtom->getId()].get(currentField); | 
|---|
| 366 | } else if (additionalAtomData.count(currentAtom->GetTrueFather()->getId())) { | 
|---|
| 367 | LOG(3, "Writing for type " << knownKeyNames[currentField] << " stuff from father: " << additionalAtomData[currentAtom->GetTrueFather()->getId()].get(currentField)); | 
|---|
| 368 | *file << additionalAtomData[currentAtom->GetTrueFather()->getId()].get(currentField); | 
|---|
| 369 | } else { | 
|---|
| 370 | LOG(3, "Writing for type " << knownKeyNames[currentField] << " AtomInfo::charge : " << currentAtom->getCharge()); | 
|---|
| 371 | *file << currentAtom->getCharge(); | 
|---|
| 372 | } | 
|---|
| 373 | } else { | 
|---|
| 374 | LOG(3, "Writing for type " << knownKeyNames[currentField] << " AtomInfo::charge : " << currentAtom->getCharge()); | 
|---|
| 375 | *file << currentAtom->getCharge(); | 
|---|
| 376 | } | 
|---|
| 377 | *file << "\t"; | 
|---|
| 378 | break; | 
|---|
| 379 | default : | 
|---|
| 380 | if (additionalAtomData.count(currentAtom->getId())) { | 
|---|
| 381 | LOG(3, "Writing for type " << knownKeyNames[currentField] << ": " << additionalAtomData[currentAtom->getId()].get(currentField)); | 
|---|
| 382 | *file << additionalAtomData[currentAtom->getId()].get(currentField); | 
|---|
| 383 | } else if (additionalAtomData.count(currentAtom->GetTrueFather()->getId())) { | 
|---|
| 384 | LOG(3, "Writing for type " << knownKeyNames[currentField] << " stuff from father: " << additionalAtomData[currentAtom->GetTrueFather()->getId()].get(currentField)); | 
|---|
| 385 | *file << additionalAtomData[currentAtom->GetTrueFather()->getId()].get(currentField); | 
|---|
| 386 | } else { | 
|---|
| 387 | LOG(3, "Writing for type " << knownKeyNames[currentField] << " the default: " << defaultAdditionalData.get(currentField)); | 
|---|
| 388 | *file << defaultAdditionalData.get(currentField); | 
|---|
| 389 | } | 
|---|
| 390 | *file << "\t"; | 
|---|
| 391 | break; | 
|---|
| 392 | } | 
|---|
| 393 | } | 
|---|
| 394 |  | 
|---|
| 395 | *file << std::endl; | 
|---|
| 396 | } | 
|---|
| 397 |  | 
|---|
| 398 | /** | 
|---|
| 399 | * Writes the neighbor information of one atom to the provided stream. | 
|---|
| 400 | * | 
|---|
| 401 | * Note that ListOfBonds of WorldTime::CurrentTime is used. | 
|---|
| 402 | * | 
|---|
| 403 | * \param stream where to write neighbor information to | 
|---|
| 404 | * \param number of neighbors | 
|---|
| 405 | * \param reference to the atom of which to take the neighbor information | 
|---|
| 406 | */ | 
|---|
| 407 | void FormatParser< tremolo >::writeNeighbors(std::ostream* file, const int numberOfNeighbors, const atom* currentAtom) { | 
|---|
| 408 | const BondList& ListOfBonds = currentAtom->getListOfBonds(); | 
|---|
| 409 | // sort bonded indices | 
|---|
| 410 | typedef std::set<atomId_t> sortedIndices; | 
|---|
| 411 | sortedIndices sortedBonds; | 
|---|
| 412 | for (BondList::const_iterator iter = ListOfBonds.begin(); | 
|---|
| 413 | iter != ListOfBonds.end(); ++iter) | 
|---|
| 414 | sortedBonds.insert(getLocalId((*iter)->GetOtherAtom(currentAtom)->getId())); | 
|---|
| 415 | // print indices | 
|---|
| 416 | sortedIndices::const_iterator currentBond = sortedBonds.begin(); | 
|---|
| 417 | for (int i = 0; i < numberOfNeighbors; i++) { | 
|---|
| 418 | *file << (currentBond != sortedBonds.end() ? (*currentBond) : 0) << "\t"; | 
|---|
| 419 | if (currentBond != sortedBonds.end()) | 
|---|
| 420 | ++currentBond; | 
|---|
| 421 | } | 
|---|
| 422 | } | 
|---|
| 423 |  | 
|---|
| 424 | /** | 
|---|
| 425 | * Stores keys from the ATOMDATA line in \a fields. | 
|---|
| 426 | * | 
|---|
| 427 | * \param line to parse the keys from | 
|---|
| 428 | * \param offset with which offset the keys begin within the line | 
|---|
| 429 | * \param fields which usedFields to use | 
|---|
| 430 | */ | 
|---|
| 431 | void FormatParser< tremolo >::parseAtomDataKeysLine( | 
|---|
| 432 | const std::string &line, | 
|---|
| 433 | const int offset, | 
|---|
| 434 | usedFields_t &fields) { | 
|---|
| 435 | std::string keyword; | 
|---|
| 436 | std::stringstream lineStream; | 
|---|
| 437 |  | 
|---|
| 438 | lineStream << line.substr(offset); | 
|---|
| 439 | while (lineStream.good()) { | 
|---|
| 440 | lineStream >> keyword; | 
|---|
| 441 | if (knownKeys[keyword.substr(0, keyword.find("="))] == TremoloKey::noKey) { | 
|---|
| 442 | // TODO: throw exception about unknown key | 
|---|
| 443 | cout << "Unknown key: " << keyword << " is not part of the tremolo format specification." << endl; | 
|---|
| 444 | break; | 
|---|
| 445 | } | 
|---|
| 446 | fields.push_back(keyword); | 
|---|
| 447 | } | 
|---|
| 448 | //LOG(1, "INFO: " << fields); | 
|---|
| 449 | } | 
|---|
| 450 |  | 
|---|
| 451 | /** Sets the properties per atom to print to .data file by parsing line from | 
|---|
| 452 | *  \a atomdata_string. | 
|---|
| 453 | * | 
|---|
| 454 | *  We just call \sa  FormatParser< tremolo >::parseAtomDataKeysLine(), however | 
|---|
| 455 | *  we clear FormatParser< tremolo >::usedFields_save. | 
|---|
| 456 | * | 
|---|
| 457 | * @param atomdata_string line to parse with space-separated values | 
|---|
| 458 | */ | 
|---|
| 459 | void FormatParser< tremolo >::setAtomData(const std::string &atomdata_string) | 
|---|
| 460 | { | 
|---|
| 461 | usedFields_save.clear(); | 
|---|
| 462 | parseAtomDataKeysLine(atomdata_string, 0, usedFields_save); | 
|---|
| 463 | } | 
|---|
| 464 |  | 
|---|
| 465 |  | 
|---|
| 466 | /** | 
|---|
| 467 | * Reads one data line of a tremolo file and interprets it according to the keys | 
|---|
| 468 | * obtained from the ATOMDATA line. | 
|---|
| 469 | * | 
|---|
| 470 | * \param line to parse as an atom | 
|---|
| 471 | * \param *newmol molecule to add atom to | 
|---|
| 472 | */ | 
|---|
| 473 | void FormatParser< tremolo >::readAtomDataLine(const std::string &line, molecule *newmol = NULL) { | 
|---|
| 474 | std::stringstream lineStream; | 
|---|
| 475 | atom* newAtom = World::getInstance().createAtom(); | 
|---|
| 476 | const atomId_t atomid = newAtom->getId(); | 
|---|
| 477 | additionalAtomData[atomid] = TremoloAtomInfoContainer(); // fill with default values | 
|---|
| 478 | TremoloAtomInfoContainer *atomInfo = &additionalAtomData[atomid]; | 
|---|
| 479 | TremoloKey::atomDataKey currentField; | 
|---|
| 480 | ConvertTo<double> toDouble; | 
|---|
| 481 | ConvertTo<int> toInt; | 
|---|
| 482 | Vector tempVector; | 
|---|
| 483 |  | 
|---|
| 484 | // setup tokenizer, splitting up white-spaced entries | 
|---|
| 485 | typedef boost::tokenizer<boost::char_separator<char> > | 
|---|
| 486 | tokenizer; | 
|---|
| 487 | boost::char_separator<char> whitespacesep(" \t"); | 
|---|
| 488 | tokenizer tokens(line, whitespacesep); | 
|---|
| 489 | ASSERT(tokens.begin() != tokens.end(), | 
|---|
| 490 | "FormatParser< tremolo >::readAtomDataLine - empty string, need at least ' '!"); | 
|---|
| 491 | tokenizer::iterator tok_iter = tokens.begin(); | 
|---|
| 492 | // then associate each token to each file | 
|---|
| 493 | for (usedFields_t::const_iterator it = usedFields_load.begin(); it < usedFields_load.end(); it++) { | 
|---|
| 494 | const std::string keyName = it->substr(0, it->find("=")); | 
|---|
| 495 | currentField = knownKeys[keyName]; | 
|---|
| 496 | const std::string word = *tok_iter; | 
|---|
| 497 | LOG(4, "INFO: Parsing key " << keyName << " with remaining data " << word); | 
|---|
| 498 | switch (currentField) { | 
|---|
| 499 | case TremoloKey::x : | 
|---|
| 500 | // for the moment, assume there are always three dimensions | 
|---|
| 501 | for (int i=0;i<NDIM;i++) { | 
|---|
| 502 | ASSERT(tok_iter != tokens.end(), "FormatParser< tremolo >::readAtomDataLine() - no value for x["+toString(i)+"]!"); | 
|---|
| 503 | LOG(4, "INFO: Parsing key " << keyName << " with next token " << *tok_iter); | 
|---|
| 504 | newAtom->set(i, toDouble(*tok_iter)); | 
|---|
| 505 | tok_iter++; | 
|---|
| 506 | } | 
|---|
| 507 | break; | 
|---|
| 508 | case TremoloKey::u : | 
|---|
| 509 | // for the moment, assume there are always three dimensions | 
|---|
| 510 | for (int i=0;i<NDIM;i++) { | 
|---|
| 511 | ASSERT(tok_iter != tokens.end(), "FormatParser< tremolo >::readAtomDataLine() - no value for u["+toString(i)+"]!"); | 
|---|
| 512 | LOG(4, "INFO: Parsing key " << keyName << " with next token " << *tok_iter); | 
|---|
| 513 | tempVector[i] = toDouble(*tok_iter); | 
|---|
| 514 | tok_iter++; | 
|---|
| 515 | } | 
|---|
| 516 | newAtom->setAtomicVelocity(tempVector); | 
|---|
| 517 | break; | 
|---|
| 518 | case TremoloKey::type : | 
|---|
| 519 | { | 
|---|
| 520 | ASSERT(tok_iter != tokens.end(), "FormatParser< tremolo >::readAtomDataLine() - no value for "+keyName+"!"); | 
|---|
| 521 | LOG(4, "INFO: Parsing key " << keyName << " with next token " << *tok_iter); | 
|---|
| 522 | std::string element; | 
|---|
| 523 | try { | 
|---|
| 524 | element = knownTypes.getType(*tok_iter); | 
|---|
| 525 | } catch(IllegalParserKeyException) { | 
|---|
| 526 | // clean up | 
|---|
| 527 | World::getInstance().destroyAtom(newAtom); | 
|---|
| 528 | // give an error | 
|---|
| 529 | ELOG(0, "TremoloParser: I do not understand the element token " << *tok_iter << "."); | 
|---|
| 530 | } | 
|---|
| 531 | // put type name into container for later use | 
|---|
| 532 | atomInfo->set(currentField, *tok_iter); | 
|---|
| 533 | LOG(4, "INFO: Parsing element " << (*tok_iter) << " as " << element << " according to KnownTypes."); | 
|---|
| 534 | tok_iter++; | 
|---|
| 535 | newAtom->setType(World::getInstance().getPeriode()->FindElement(element)); | 
|---|
| 536 | ASSERT(newAtom->getType(), "Type was not set for this atom"); | 
|---|
| 537 | break; | 
|---|
| 538 | } | 
|---|
| 539 | case TremoloKey::Id : | 
|---|
| 540 | ASSERT(tok_iter != tokens.end(), "FormatParser< tremolo >::readAtomDataLine() - no value for "+keyName+"!"); | 
|---|
| 541 | LOG(4, "INFO: Parsing key " << keyName << " with next token " << *tok_iter); | 
|---|
| 542 | associateLocaltoGlobalId(toInt(*tok_iter), atomid); | 
|---|
| 543 | tok_iter++; | 
|---|
| 544 | break; | 
|---|
| 545 | case TremoloKey::neighbors : | 
|---|
| 546 | for (int i=0;i<atoi(it->substr(it->find("=") + 1, 1).c_str());i++) { | 
|---|
| 547 | ASSERT(tok_iter != tokens.end(), "FormatParser< tremolo >::readAtomDataLine() - no value for "+keyName+"!"); | 
|---|
| 548 | LOG(4, "INFO: Parsing key " << keyName << " with next token " << *tok_iter); | 
|---|
| 549 | lineStream << *tok_iter << "\t"; | 
|---|
| 550 | tok_iter++; | 
|---|
| 551 | } | 
|---|
| 552 | readNeighbors(&lineStream, | 
|---|
| 553 | atoi(it->substr(it->find("=") + 1, 1).c_str()), atomid); | 
|---|
| 554 | break; | 
|---|
| 555 | case TremoloKey::charge : | 
|---|
| 556 | ASSERT(tok_iter != tokens.end(), "FormatParser< tremolo >::readAtomDataLine() - no value for "+keyName+"!"); | 
|---|
| 557 | LOG(4, "INFO: Parsing key " << keyName << " with next token " << *tok_iter); | 
|---|
| 558 | atomInfo->set(currentField, *tok_iter); | 
|---|
| 559 | newAtom->setCharge(boost::lexical_cast<double>(*tok_iter)); | 
|---|
| 560 | tok_iter++; | 
|---|
| 561 | break; | 
|---|
| 562 | default : | 
|---|
| 563 | ASSERT(tok_iter != tokens.end(), "FormatParser< tremolo >::readAtomDataLine() - no value for "+keyName+"!"); | 
|---|
| 564 | LOG(4, "INFO: Parsing key " << keyName << " with next token " << *tok_iter); | 
|---|
| 565 | atomInfo->set(currentField, *tok_iter); | 
|---|
| 566 | tok_iter++; | 
|---|
| 567 | break; | 
|---|
| 568 | } | 
|---|
| 569 | } | 
|---|
| 570 | LOG(3, "INFO: Parsed atom " << atomid << "."); | 
|---|
| 571 | if (newmol != NULL) | 
|---|
| 572 | newmol->AddAtom(newAtom); | 
|---|
| 573 | } | 
|---|
| 574 |  | 
|---|
| 575 | bool FormatParser< tremolo >::saveAtomsInExttypes(std::ostream &output, const std::vector<atom*> &atoms, const int id) const | 
|---|
| 576 | { | 
|---|
| 577 | bool status = true; | 
|---|
| 578 | // parse the file | 
|---|
| 579 | for (std::vector<atom *>::const_iterator iter = atoms.begin(); | 
|---|
| 580 | iter != atoms.end(); ++iter) { | 
|---|
| 581 | const int atomicid = getLocalId((*iter)->getId()); | 
|---|
| 582 | if (atomicid == -1) | 
|---|
| 583 | status = false; | 
|---|
| 584 | output << atomicid << "\t" << id << std::endl; | 
|---|
| 585 | } | 
|---|
| 586 |  | 
|---|
| 587 | return status; | 
|---|
| 588 | } | 
|---|
| 589 |  | 
|---|
| 590 | /** | 
|---|
| 591 | * Reads neighbor information for one atom from the input. | 
|---|
| 592 | * | 
|---|
| 593 | * \param line stream where to read the information from | 
|---|
| 594 | * \param numberOfNeighbors number of neighbors to read | 
|---|
| 595 | * \param atomid world id of the atom the information belongs to | 
|---|
| 596 | */ | 
|---|
| 597 | void FormatParser< tremolo >::readNeighbors(std::stringstream* line, const int numberOfNeighbors, const int atomId) { | 
|---|
| 598 | int neighborId = 0; | 
|---|
| 599 | for (int i = 0; i < numberOfNeighbors; i++) { | 
|---|
| 600 | *line >> neighborId; | 
|---|
| 601 | // 0 is used to fill empty neighbor positions in the tremolo file. | 
|---|
| 602 | if (neighborId > 0) { | 
|---|
| 603 | LOG(4, "INFO: Atom with global id " << atomId | 
|---|
| 604 | << " has neighbour with serial " << neighborId); | 
|---|
| 605 | additionalAtomData[atomId].neighbors.push_back(neighborId); | 
|---|
| 606 | } | 
|---|
| 607 | } | 
|---|
| 608 | } | 
|---|
| 609 |  | 
|---|
| 610 | /** | 
|---|
| 611 | * Checks whether the provided name is within \a fields. | 
|---|
| 612 | * | 
|---|
| 613 | * \param fields which usedFields to use | 
|---|
| 614 | * \param fieldName name to check | 
|---|
| 615 | * \return true if the field name is used | 
|---|
| 616 | */ | 
|---|
| 617 | bool FormatParser< tremolo >::isUsedField(const usedFields_t &fields, const std::string &fieldName) const | 
|---|
| 618 | { | 
|---|
| 619 | bool fieldNameExists = false; | 
|---|
| 620 | for (usedFields_t::const_iterator usedField = fields.begin(); | 
|---|
| 621 | usedField != fields.end(); usedField++) { | 
|---|
| 622 | if (usedField->substr(0, usedField->find("=")) == fieldName) | 
|---|
| 623 | fieldNameExists = true; | 
|---|
| 624 | } | 
|---|
| 625 |  | 
|---|
| 626 | return fieldNameExists; | 
|---|
| 627 | } | 
|---|
| 628 |  | 
|---|
| 629 |  | 
|---|
| 630 | /** | 
|---|
| 631 | * Adds the collected neighbor information to the atoms in the world. The atoms | 
|---|
| 632 | * are found by their current ID and mapped to the corresponding atoms with the | 
|---|
| 633 | * Id found in the parsed file. | 
|---|
| 634 | * | 
|---|
| 635 | * @param atoms vector with all newly added (global) atomic ids | 
|---|
| 636 | */ | 
|---|
| 637 | void FormatParser< tremolo >::processNeighborInformation(const std::vector<atomId_t> &atoms) { | 
|---|
| 638 | if (!isUsedField(usedFields_load, "neighbors")) { | 
|---|
| 639 | return; | 
|---|
| 640 | } | 
|---|
| 641 |  | 
|---|
| 642 | for (std::vector<atomId_t>::const_iterator iter = atoms.begin(); iter != atoms.end(); ++iter) { | 
|---|
| 643 | ASSERT(additionalAtomData.count(*iter) != 0, | 
|---|
| 644 | "FormatParser< tremolo >::processNeighborInformation() - global id " | 
|---|
| 645 | +toString(*iter)+" unknown in additionalAtomData."); | 
|---|
| 646 | TremoloAtomInfoContainer ¤tInfo = additionalAtomData[*iter]; | 
|---|
| 647 | ASSERT (!currentInfo.neighbors_processed, | 
|---|
| 648 | "FormatParser< tremolo >::processNeighborInformation() - neighbors of new atom " | 
|---|
| 649 | +toString(*iter)+" are already processed."); | 
|---|
| 650 | for(std::vector<int>::const_iterator neighbor = currentInfo.neighbors.begin(); | 
|---|
| 651 | neighbor != currentInfo.neighbors.end(); neighbor++ | 
|---|
| 652 | ) { | 
|---|
| 653 | LOG(3, "INFO: Creating bond between (" | 
|---|
| 654 | << *iter | 
|---|
| 655 | << ") and (" | 
|---|
| 656 | << getGlobalId(*neighbor) << "|" << *neighbor << ")"); | 
|---|
| 657 | ASSERT(getGlobalId(*neighbor) != -1, | 
|---|
| 658 | "FormatParser< tremolo >::processNeighborInformation() - global id to local id " | 
|---|
| 659 | +toString(*neighbor)+" is unknown."); | 
|---|
| 660 | World::getInstance().getAtom(AtomById(*iter)) | 
|---|
| 661 | ->addBond(WorldTime::getTime(), World::getInstance().getAtom(AtomById(getGlobalId(*neighbor)))); | 
|---|
| 662 | } | 
|---|
| 663 | currentInfo.neighbors_processed = true; | 
|---|
| 664 | } | 
|---|
| 665 | } | 
|---|
| 666 |  | 
|---|
| 667 | /** | 
|---|
| 668 | * Replaces atom IDs read from the file by the corresponding world IDs. All IDs | 
|---|
| 669 | * IDs of the input string will be replaced; expected separating characters are | 
|---|
| 670 | * "-" and ",". | 
|---|
| 671 | * | 
|---|
| 672 | * \param string in which atom IDs should be adapted | 
|---|
| 673 | * | 
|---|
| 674 | * \return input string with modified atom IDs | 
|---|
| 675 | */ | 
|---|
| 676 | std::string FormatParser< tremolo >::adaptIdDependentDataString(std::string data) { | 
|---|
| 677 | // there might be no IDs | 
|---|
| 678 | if (data == "-") { | 
|---|
| 679 | return "-"; | 
|---|
| 680 | } | 
|---|
| 681 |  | 
|---|
| 682 | char separator; | 
|---|
| 683 | int id; | 
|---|
| 684 | std::stringstream line, result; | 
|---|
| 685 | line << data; | 
|---|
| 686 |  | 
|---|
| 687 | line >> id; | 
|---|
| 688 | result << getGlobalId(id); | 
|---|
| 689 | while (line.good()) { | 
|---|
| 690 | line >> separator >> id; | 
|---|
| 691 | result << separator << getGlobalId(id); | 
|---|
| 692 | } | 
|---|
| 693 |  | 
|---|
| 694 | return result.str(); | 
|---|
| 695 | } | 
|---|
| 696 |  | 
|---|
| 697 | /** | 
|---|
| 698 | * Corrects the atom IDs in each imprData entry to the corresponding world IDs | 
|---|
| 699 | * as they might differ from the originally read IDs. | 
|---|
| 700 | */ | 
|---|
| 701 | void FormatParser< tremolo >::adaptImprData() { | 
|---|
| 702 | if (!isUsedField(usedFields_load, "imprData")) { | 
|---|
| 703 | return; | 
|---|
| 704 | } | 
|---|
| 705 |  | 
|---|
| 706 | for(std::map<const atomId_t, TremoloAtomInfoContainer>::iterator currentInfo = additionalAtomData.begin(); | 
|---|
| 707 | currentInfo != additionalAtomData.end(); currentInfo++ | 
|---|
| 708 | ) { | 
|---|
| 709 | currentInfo->second.imprData = adaptIdDependentDataString(currentInfo->second.imprData); | 
|---|
| 710 | } | 
|---|
| 711 | } | 
|---|
| 712 |  | 
|---|
| 713 | /** | 
|---|
| 714 | * Corrects the atom IDs in each torsion entry to the corresponding world IDs | 
|---|
| 715 | * as they might differ from the originally read IDs. | 
|---|
| 716 | */ | 
|---|
| 717 | void FormatParser< tremolo >::adaptTorsion() { | 
|---|
| 718 | if (!isUsedField(usedFields_load, "torsion")) { | 
|---|
| 719 | return; | 
|---|
| 720 | } | 
|---|
| 721 |  | 
|---|
| 722 | for(std::map<const atomId_t, TremoloAtomInfoContainer>::iterator currentInfo = additionalAtomData.begin(); | 
|---|
| 723 | currentInfo != additionalAtomData.end(); currentInfo++ | 
|---|
| 724 | ) { | 
|---|
| 725 | currentInfo->second.torsion = adaptIdDependentDataString(currentInfo->second.torsion); | 
|---|
| 726 | } | 
|---|
| 727 | } | 
|---|
| 728 |  | 
|---|