/*
* 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++;
}
}
}
}