/* * Project: MoleCuilder * Description: creates and alters molecular systems * Copyright (C) 2011 University of Bonn. All rights reserved. * * * This file is part of MoleCuilder. * * MoleCuilder is free software: you can redistribute it and/or modify * it under the terms of the GNU General Public License as published by * the Free Software Foundation, either version 2 of the License, or * (at your option) any later version. * * MoleCuilder is distributed in the hope that it will be useful, * but WITHOUT ANY WARRANTY; without even the implied warranty of * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the * GNU General Public License for more details. * * You should have received a copy of the GNU General Public License * along with MoleCuilder. If not, see . */ /* * ExportGraph.cpp * * Created on: 08.03.2012 * Author: heber */ // include config.h #ifdef HAVE_CONFIG_H #include #endif #include "CodePatterns/MemDebug.hpp" #include "ExportGraph.hpp" #include "CodePatterns/Info.hpp" #include "CodePatterns/Log.hpp" #include "Bond/bond.hpp" #include "Element/element.hpp" #include "Fragmentation/Graph.hpp" #include "Fragmentation/KeySet.hpp" #include "Fragmentation/SortIndex.hpp" #include "Graph/ListOfLocalAtoms.hpp" #include "molecule.hpp" #include "MoleculeListClass.hpp" #include "World.hpp" /** Constructor for class ExportGraph. * * @param _graph */ ExportGraph::ExportGraph( const Graph &_graph, const enum HydrogenTreatment _treatment, const enum HydrogenSaturation _saturation) : TotalGraph(_graph), CurrentKeySet(TotalGraph.begin()), BondFragments(World::getPointer()), treatment(_treatment), saturation(_saturation) { } /** Destructor for class ExportGraph. * */ ExportGraph::~ExportGraph() { // remove all create molecules again from the World including their atoms for (MoleculeList::iterator iter = BondFragments.ListOfMolecules.begin(); !BondFragments.ListOfMolecules.empty(); iter = BondFragments.ListOfMolecules.begin()) { // remove copied atoms and molecule again molecule *mol = *iter; mol->removeAtomsinMolecule(); World::getInstance().destroyMolecule(mol); BondFragments.ListOfMolecules.erase(iter); } } void ExportGraph::operator()() { if (BondFragments.ListOfMolecules.size() == 0) prepareMolecule(); } void ExportGraph::reset() { CurrentKeySet = TotalGraph.begin(); } ExportGraph::SaturatedFragment_ptr ExportGraph::getNextFragment() { // if a fragment is still leased, return zero ptr. if (!KeySetsInUse.empty()) { ELOG(1, "Leasing KeySet while old one is not returned."); return SaturatedFragment_ptr(); } // else return current fragment or indicate end if (CurrentKeySet != TotalGraph.end()) { const KeySet &set = (CurrentKeySet++)->first; return leaseFragment(set); } else { return leaseFragment(EmptySet); } } ExportGraph::SaturatedFragment_ptr ExportGraph::leaseFragment(const KeySet &_set) { // create the saturation which adds itself to KeySetsInUse SaturatedFragment_ptr _ptr( new SaturatedFragment( _set, KeySetsInUse, hydrogens, treatment, saturation) ); // and return return _ptr; } void ExportGraph::releaseFragment(SaturatedFragment_ptr &_ptr) { ASSERT(_ptr != NULL, "ExportGraph::releaseFragment() - pointer is NULL."); SaturatedFragment::KeySetsInUse_t::iterator iter = KeySetsInUse.find(_ptr->getKeySet()); if (iter == KeySetsInUse.end()) { ASSERT(0, "ExportGraph::releaseFragment() - returning unknown set " +toString(_ptr->getKeySet())+"."); return; } else { // release instance which removes itself in KeySetsInUse _ptr.reset(); } } /** Internal helper to create from each keyset a molecule * */ void ExportGraph::prepareMolecule() { size_t count = 0; for(Graph::const_iterator runner = TotalGraph.begin(); runner != TotalGraph.end(); runner++) { KeySet test = (*runner).first; LOG(2, "DEBUG: Fragment No." << (*runner).second.first << " with TEFactor " << (*runner).second.second << "."); BondFragments.insert(StoreFragmentFromKeySet(test, World::getInstance().getConfig())); ++count; } LOG(1, "INFO: " << count << "/" << BondFragments.ListOfMolecules.size() << " fragments generated from the keysets."); } /** Stores a fragment from \a KeySet into \a molecule. * First creates the minimal set of atoms from the KeySet, then creates the bond structure from the complete * molecule and adds missing hydrogen where bonds were cut. * \param &Leaflet pointer to KeySet structure * \param IsAngstroem whether we have Ansgtroem or bohrradius * \return pointer to constructed molecule */ molecule * ExportGraph::StoreFragmentFromKeySet(KeySet &Leaflet, bool IsAngstroem) { Info info(__func__); ListOfLocalAtoms_t SonList; molecule *Leaf = World::getInstance().createMolecule(); StoreFragmentFromKeySet_Init(Leaf, Leaflet, SonList); // create the bonds between all: Make it an induced subgraph and add hydrogen // LOG(2, "Creating bonds from father graph (i.e. induced subgraph creation)."); CreateInducedSubgraphOfFragment(Leaf, SonList, IsAngstroem); //Leaflet->Leaf->ScanForPeriodicCorrection(out); return Leaf; } /** Initializes some value for putting fragment of \a *mol into \a *Leaf. * \param *Leaf fragment molecule * \param &Leaflet pointer to KeySet structure * \param SonList calloc'd list which atom of \a *Leaf is a son of which atom in \a *mol * \return number of atoms in fragment */ int ExportGraph::StoreFragmentFromKeySet_Init(molecule *Leaf, KeySet &Leaflet, ListOfLocalAtoms_t &SonList) { atom *FatherOfRunner = NULL; // first create the minimal set of atoms from the KeySet World &world = World::getInstance(); int size = 0; for(KeySet::const_iterator runner = Leaflet.begin(); runner != Leaflet.end(); runner++) { FatherOfRunner = world.getAtom(AtomById(*runner)); // find the id SonList.insert( std::make_pair(FatherOfRunner->getNr(), Leaf->AddCopyAtom(FatherOfRunner) ) ); size++; } return size; } /** Creates an induced subgraph out of a fragmental key set, adding bonds and hydrogens (if treated specially). * \param *Leaf fragment molecule * \param IsAngstroem whether we have Ansgtroem or bohrradius * \param SonList list which atom of \a *Leaf is another atom's son */ void ExportGraph::CreateInducedSubgraphOfFragment(molecule *Leaf, ListOfLocalAtoms_t &SonList, bool IsAngstroem) { bool LonelyFlag = false; atom *OtherFather = NULL; atom *FatherOfRunner = NULL; // we increment the iter just before skipping the hydrogen // as we use AddBond, we cannot have a const_iterator here for (molecule::iterator iter = Leaf->begin(); iter != Leaf->end();) { LonelyFlag = true; FatherOfRunner = (*iter)->father; ASSERT(FatherOfRunner,"Atom without father found"); if (SonList.find(FatherOfRunner->getNr()) != SonList.end()) { // check if this, our father, is present in list // create all bonds const BondList& ListOfBonds = FatherOfRunner->getListOfBonds(); for (BondList::const_iterator BondRunner = ListOfBonds.begin(); BondRunner != ListOfBonds.end(); ++BondRunner) { OtherFather = (*BondRunner)->GetOtherAtom(FatherOfRunner); if (SonList.find(OtherFather->getNr()) != SonList.end()) { // LOG(2, "INFO: Father " << *FatherOfRunner << " of son " << *SonList[FatherOfRunner->getNr()] // << " is bound to " << *OtherFather << ", whose son is " // << *SonList[OtherFather->getNr()] << "."); if (OtherFather->getNr() > FatherOfRunner->getNr()) { // add bond (Nr check is for adding only one of both variants: ab, ba) std::stringstream output; // output << "ACCEPT: Adding Bond: " output << Leaf->AddBond((*iter), SonList[OtherFather->getNr()], (*BondRunner)->BondDegree); // LOG(3, output.str()); //NumBonds[(*iter)->getNr()]++; } else { // LOG(3, "REJECY: Not adding bond, labels in wrong order."); } LonelyFlag = false; } else { // LOG(2, "INFO: Father " << *FatherOfRunner << " of son " << *SonList[FatherOfRunner->getNr()] // << " is bound to " << *OtherFather << ", who has no son in this fragment molecule."); if (saturation == DoSaturate) { // LOG(3, "ACCEPT: Adding Hydrogen to " << (*iter)->Name << " and a bond in between."); if (!Leaf->AddHydrogenReplacementAtom((*BondRunner), (*iter), FatherOfRunner, OtherFather, IsAngstroem)) exit(1); } else if ((treatment == ExcludeHydrogen) && (OtherFather->getElementNo() == (atomicNumber_t)1)) { // just copy the atom if it's a hydrogen atom * const OtherWalker = Leaf->AddCopyAtom(OtherFather); Leaf->AddBond((*iter), OtherWalker, (*BondRunner)->BondDegree); } //NumBonds[(*iter)->getNr()] += Binder->BondDegree; } } } else { ELOG(1, "Son " << (*iter)->getName() << " has father " << FatherOfRunner->getName() << " but its entry in SonList is " << SonList[FatherOfRunner->getNr()] << "!"); } if ((LonelyFlag) && (Leaf->getAtomCount() > 1)) { LOG(0, **iter << "has got bonds only to hydrogens!"); } ++iter; if (saturation == DoSaturate) { while ((iter != Leaf->end()) && ((*iter)->getType()->getAtomicNumber() == 1)){ // skip added hydrogen iter++; } } } }