source: src/Parser/TremoloParser.cpp@ ca331c

Action_Thermostats Add_AtomRandomPerturbation Add_FitFragmentPartialChargesAction Add_RotateAroundBondAction Add_SelectAtomByNameAction Added_ParseSaveFragmentResults AddingActions_SaveParseParticleParameters Adding_Graph_to_ChangeBondActions Adding_MD_integration_tests Adding_ParticleName_to_Atom Adding_StructOpt_integration_tests AtomFragments Automaking_mpqc_open AutomationFragmentation_failures Candidate_v1.5.4 Candidate_v1.6.0 Candidate_v1.6.1 ChangeBugEmailaddress ChangingTestPorts ChemicalSpaceEvaluator CombiningParticlePotentialParsing Combining_Subpackages Debian_Package_split Debian_package_split_molecuildergui_only Disabling_MemDebug Docu_Python_wait EmpiricalPotential_contain_HomologyGraph EmpiricalPotential_contain_HomologyGraph_documentation Enable_parallel_make_install Enhance_userguide Enhanced_StructuralOptimization Enhanced_StructuralOptimization_continued Example_ManyWaysToTranslateAtom Exclude_Hydrogens_annealWithBondGraph FitPartialCharges_GlobalError Fix_BoundInBox_CenterInBox_MoleculeActions Fix_ChargeSampling_PBC Fix_ChronosMutex Fix_FitPartialCharges Fix_FitPotential_needs_atomicnumbers Fix_ForceAnnealing Fix_IndependentFragmentGrids Fix_ParseParticles Fix_ParseParticles_split_forward_backward_Actions Fix_PopActions Fix_QtFragmentList_sorted_selection Fix_Restrictedkeyset_FragmentMolecule Fix_StatusMsg Fix_StepWorldTime_single_argument Fix_Verbose_Codepatterns Fix_fitting_potentials Fixes ForceAnnealing_goodresults ForceAnnealing_oldresults ForceAnnealing_tocheck ForceAnnealing_with_BondGraph ForceAnnealing_with_BondGraph_continued ForceAnnealing_with_BondGraph_continued_betteresults ForceAnnealing_with_BondGraph_contraction-expansion FragmentAction_writes_AtomFragments FragmentMolecule_checks_bonddegrees GeometryObjects Gui_Fixes Gui_displays_atomic_force_velocity ImplicitCharges IndependentFragmentGrids IndependentFragmentGrids_IndividualZeroInstances IndependentFragmentGrids_IntegrationTest IndependentFragmentGrids_Sole_NN_Calculation JobMarket_RobustOnKillsSegFaults JobMarket_StableWorkerPool JobMarket_unresolvable_hostname_fix MoreRobust_FragmentAutomation ODR_violation_mpqc_open PartialCharges_OrthogonalSummation PdbParser_setsAtomName PythonUI_with_named_parameters QtGui_reactivate_TimeChanged_changes Recreated_GuiChecks Rewrite_FitPartialCharges RotateToPrincipalAxisSystem_UndoRedo SaturateAtoms_findBestMatching SaturateAtoms_singleDegree StoppableMakroAction Subpackage_CodePatterns Subpackage_JobMarket Subpackage_LinearAlgebra Subpackage_levmar Subpackage_mpqc_open Subpackage_vmg Switchable_LogView ThirdParty_MPQC_rebuilt_buildsystem TrajectoryDependenant_MaxOrder TremoloParser_IncreasedPrecision TremoloParser_MultipleTimesteps TremoloParser_setsAtomName Ubuntu_1604_changes stable
Last change on this file since ca331c was ca331c, checked in by Frederik Heber <heber@…>, 13 years ago

Added Undo/Redo capability to SetTremoloAtomdataAction.

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