| [dadc74] | 1 | /*
 | 
|---|
 | 2 |  * Project: MoleCuilder
 | 
|---|
 | 3 |  * Description: creates and alters molecular systems
 | 
|---|
| [0aa122] | 4 |  * Copyright (C)  2010-2012 University of Bonn. All rights reserved.
 | 
|---|
| [94d5ac6] | 5 |  * 
 | 
|---|
 | 6 |  *
 | 
|---|
 | 7 |  *   This file is part of MoleCuilder.
 | 
|---|
 | 8 |  *
 | 
|---|
 | 9 |  *    MoleCuilder is free software: you can redistribute it and/or modify
 | 
|---|
 | 10 |  *    it under the terms of the GNU General Public License as published by
 | 
|---|
 | 11 |  *    the Free Software Foundation, either version 2 of the License, or
 | 
|---|
 | 12 |  *    (at your option) any later version.
 | 
|---|
 | 13 |  *
 | 
|---|
 | 14 |  *    MoleCuilder is distributed in the hope that it will be useful,
 | 
|---|
 | 15 |  *    but WITHOUT ANY WARRANTY; without even the implied warranty of
 | 
|---|
 | 16 |  *    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 | 
|---|
 | 17 |  *    GNU General Public License for more details.
 | 
|---|
 | 18 |  *
 | 
|---|
 | 19 |  *    You should have received a copy of the GNU General Public License
 | 
|---|
 | 20 |  *    along with MoleCuilder.  If not, see <http://www.gnu.org/licenses/>.
 | 
|---|
| [dadc74] | 21 |  */
 | 
|---|
 | 22 | 
 | 
|---|
 | 23 | /*
 | 
|---|
 | 24 |  * Graph.cpp
 | 
|---|
 | 25 |  *
 | 
|---|
 | 26 |  *  Created on: Oct 20, 2011
 | 
|---|
 | 27 |  *      Author: heber
 | 
|---|
 | 28 |  */
 | 
|---|
 | 29 | 
 | 
|---|
 | 30 | // include config.h
 | 
|---|
 | 31 | #ifdef HAVE_CONFIG_H
 | 
|---|
 | 32 | #include <config.h>
 | 
|---|
 | 33 | #endif
 | 
|---|
 | 34 | 
 | 
|---|
 | 35 | #include "CodePatterns/MemDebug.hpp"
 | 
|---|
 | 36 | 
 | 
|---|
| [75363b] | 37 | #include <fstream>
 | 
|---|
 | 38 | #include <iostream>
 | 
|---|
 | 39 | #include <sstream>
 | 
|---|
 | 40 | #include <string>
 | 
|---|
 | 41 | 
 | 
|---|
| [dadc74] | 42 | #include "Graph.hpp"
 | 
|---|
 | 43 | 
 | 
|---|
 | 44 | #include "CodePatterns/Log.hpp"
 | 
|---|
 | 45 | 
 | 
|---|
| [bfbd4a] | 46 | #include "Atom/atom.hpp"
 | 
|---|
 | 47 | #include "Descriptors/AtomIdDescriptor.hpp"
 | 
|---|
| [730d7a] | 48 | #include "Fragmentation/AdaptivityMap.hpp"
 | 
|---|
| [75363b] | 49 | #include "Helpers/defs.hpp"
 | 
|---|
 | 50 | #include "Helpers/helpers.hpp"
 | 
|---|
| [bfbd4a] | 51 | #include "World.hpp"
 | 
|---|
| [75363b] | 52 | 
 | 
|---|
| [dadc74] | 53 | /** Constructor for class Graph.
 | 
|---|
 | 54 |  *
 | 
|---|
 | 55 |  */
 | 
|---|
 | 56 | Graph::Graph()
 | 
|---|
 | 57 | {}
 | 
|---|
 | 58 | 
 | 
|---|
 | 59 | /** Destructor for class Graph.
 | 
|---|
 | 60 |  *
 | 
|---|
 | 61 |  */
 | 
|---|
 | 62 | Graph::~Graph()
 | 
|---|
 | 63 | {}
 | 
|---|
 | 64 | 
 | 
|---|
 | 65 | /** Inserts each KeySet in \a graph into \a this.
 | 
|---|
 | 66 |  * \param graph1 graph whose KeySet are inserted into \this graph
 | 
|---|
 | 67 |  * \param *counter keyset counter that gets increased
 | 
|---|
 | 68 |  */
 | 
|---|
| [b4f72c] | 69 | void Graph::InsertGraph(const Graph &graph, int &counter)
 | 
|---|
| [dadc74] | 70 | {
 | 
|---|
 | 71 |   GraphTestPair testGraphInsert;
 | 
|---|
 | 72 | 
 | 
|---|
| [b4f72c] | 73 |   for(Graph::const_iterator runner = graph.begin(); runner != graph.end(); runner++) {
 | 
|---|
 | 74 |     testGraphInsert = insert(GraphPair ((*runner).first,pair<int,double>(++counter,((*runner).second).second)));  // store fragment number and current factor
 | 
|---|
| [dadc74] | 75 |     if (testGraphInsert.second) {
 | 
|---|
| [b4f72c] | 76 |       LOG(2, "INFO: KeySet " << counter-1 << " successfully inserted.");
 | 
|---|
| [dadc74] | 77 |     } else {
 | 
|---|
| [b4f72c] | 78 |       LOG(2, "INFO: KeySet " << counter-1 << " failed to insert, present fragment is " << ((*(testGraphInsert.first)).second).first);
 | 
|---|
| [dadc74] | 79 |       ((*(testGraphInsert.first)).second).second += (*runner).second.second;
 | 
|---|
 | 80 |       LOG(2, "INFO: New factor is " << (*(testGraphInsert.first)).second.second << ".");
 | 
|---|
 | 81 |     }
 | 
|---|
 | 82 |   }
 | 
|---|
 | 83 | };
 | 
|---|
 | 84 | 
 | 
|---|
| [75363b] | 85 | /** Parses the KeySet file and fills \a this from the known molecule structure.
 | 
|---|
 | 86 |  * Does two-pass scanning:
 | 
|---|
 | 87 |  * -# Scans the keyset file and initialises a temporary graph
 | 
|---|
 | 88 |  * -# Scans TEFactors file and sets the TEFactor of each key set in the temporary graph accordingly
 | 
|---|
 | 89 |  * Finally, the temporary graph is inserted into the given \a FragmentList for return.
 | 
|---|
 | 90 |  * \param &path path to file
 | 
|---|
 | 91 |  * \return true - parsing successfully, false - failure on parsing (FragmentList will be NULL)
 | 
|---|
 | 92 |  */
 | 
|---|
| [e71325] | 93 | bool Graph::ParseKeySetFile(const std::string &path)
 | 
|---|
| [75363b] | 94 | {
 | 
|---|
 | 95 |   bool status = true;
 | 
|---|
 | 96 |   std::ifstream InputFile;
 | 
|---|
 | 97 |   std::stringstream line;
 | 
|---|
 | 98 |   GraphTestPair testGraphInsert;
 | 
|---|
 | 99 |   int NumberOfFragments = 0;
 | 
|---|
 | 100 |   std::string filename;
 | 
|---|
 | 101 | 
 | 
|---|
 | 102 |   // 1st pass: open file and read
 | 
|---|
 | 103 |   LOG(1, "INFO: Parsing the KeySet file ... ");
 | 
|---|
 | 104 |   filename = path + KEYSETFILE;
 | 
|---|
 | 105 |   InputFile.open(filename.c_str());
 | 
|---|
 | 106 |   if (InputFile.good()) {
 | 
|---|
 | 107 |     // each line represents a new fragment
 | 
|---|
 | 108 |     char buffer[MAXSTRINGSIZE];
 | 
|---|
 | 109 |     // 1. parse keysets and insert into temp. graph
 | 
|---|
 | 110 |     while (!InputFile.eof()) {
 | 
|---|
 | 111 |       InputFile.getline(buffer, MAXSTRINGSIZE);
 | 
|---|
 | 112 |       KeySet CurrentSet;
 | 
|---|
 | 113 |       if ((strlen(buffer) > 0) && (CurrentSet.ScanBufferIntoKeySet(buffer))) {  // if at least one valid atom was added, write config
 | 
|---|
 | 114 |         testGraphInsert = insert(GraphPair (CurrentSet,pair<int,double>(NumberOfFragments++,1)));  // store fragment number and current factor
 | 
|---|
 | 115 |         if (!testGraphInsert.second) {
 | 
|---|
 | 116 |           ELOG(0, "KeySet file must be corrupt as there are two equal key sets therein!");
 | 
|---|
 | 117 |           performCriticalExit();
 | 
|---|
 | 118 |         }
 | 
|---|
 | 119 |       }
 | 
|---|
 | 120 |     }
 | 
|---|
 | 121 |     // 2. Free and done
 | 
|---|
 | 122 |     InputFile.close();
 | 
|---|
 | 123 |     InputFile.clear();
 | 
|---|
 | 124 |     LOG(1, "INFO: ... done.");
 | 
|---|
 | 125 |   } else {
 | 
|---|
 | 126 |     ELOG(1, "File " << filename << " not found.");
 | 
|---|
 | 127 |     status = false;
 | 
|---|
 | 128 |   }
 | 
|---|
 | 129 | 
 | 
|---|
 | 130 |   return status;
 | 
|---|
 | 131 | };
 | 
|---|
 | 132 | 
 | 
|---|
 | 133 | /** Stores key sets to file.
 | 
|---|
 | 134 |  * \param &path path to file
 | 
|---|
 | 135 |  * \return true - file written successfully, false - writing failed
 | 
|---|
 | 136 |  */
 | 
|---|
| [e71325] | 137 | bool Graph::StoreKeySetFile(const std::string &path) const
 | 
|---|
| [75363b] | 138 | {
 | 
|---|
 | 139 |   bool status =  true;
 | 
|---|
 | 140 |   std::string line = path + KEYSETFILE;
 | 
|---|
 | 141 |   std::ofstream output(line.c_str());
 | 
|---|
 | 142 | 
 | 
|---|
 | 143 |   // open KeySet file
 | 
|---|
 | 144 |   LOG(1, "INFO: Saving key sets of the total graph ... ");
 | 
|---|
 | 145 |   if(output.good()) {
 | 
|---|
| [ca8bea] | 146 |     for(Graph::const_iterator runner = begin(); runner != end(); runner++) {
 | 
|---|
 | 147 |       for (KeySet::const_iterator sprinter = (*runner).first.begin();sprinter != (*runner).first.end(); sprinter++) {
 | 
|---|
| [75363b] | 148 |         if (sprinter != (*runner).first.begin())
 | 
|---|
 | 149 |           output << "\t";
 | 
|---|
 | 150 |         output << *sprinter;
 | 
|---|
 | 151 |       }
 | 
|---|
 | 152 |       output << std::endl;
 | 
|---|
 | 153 |     }
 | 
|---|
 | 154 |     LOG(1, "INFO: done.");
 | 
|---|
 | 155 |   } else {
 | 
|---|
 | 156 |     ELOG(0, "Unable to open " << line << " for writing keysets!");
 | 
|---|
 | 157 |     performCriticalExit();
 | 
|---|
 | 158 |     status = false;
 | 
|---|
 | 159 |   }
 | 
|---|
 | 160 |   output.close();
 | 
|---|
 | 161 | 
 | 
|---|
 | 162 |   return status;
 | 
|---|
 | 163 | };
 | 
|---|
 | 164 | 
 | 
|---|
 | 165 | /** Parses the TE factors file and fills \a *FragmentList from the known molecule structure.
 | 
|---|
 | 166 |  * -# Scans TEFactors file and sets the TEFactor of each key set in the temporary graph accordingly
 | 
|---|
 | 167 |  * \param *path path to file
 | 
|---|
 | 168 |  * \return true - parsing successfully, false - failure on parsing
 | 
|---|
 | 169 |  */
 | 
|---|
 | 170 | bool Graph::ParseTEFactorsFile(char *path)
 | 
|---|
 | 171 | {
 | 
|---|
 | 172 |   bool status = true;
 | 
|---|
 | 173 |   std::ifstream InputFile;
 | 
|---|
 | 174 |   std::stringstream line;
 | 
|---|
 | 175 |   GraphTestPair testGraphInsert;
 | 
|---|
 | 176 |   int NumberOfFragments = 0;
 | 
|---|
 | 177 |   double TEFactor;
 | 
|---|
 | 178 |   char filename[MAXSTRINGSIZE];
 | 
|---|
 | 179 | 
 | 
|---|
 | 180 |   // 2nd pass: open TEFactors file and read
 | 
|---|
 | 181 |   LOG(1, "INFO: Parsing the TEFactors file ... ");
 | 
|---|
 | 182 |   sprintf(filename, "%s/%s%s", path, FRAGMENTPREFIX, TEFACTORSFILE);
 | 
|---|
 | 183 |   InputFile.open(filename);
 | 
|---|
 | 184 |   if (InputFile != NULL) {
 | 
|---|
 | 185 |     // 3. add found TEFactors to each keyset
 | 
|---|
 | 186 |     NumberOfFragments = 0;
 | 
|---|
 | 187 |     for(Graph::iterator runner = begin();runner != end(); runner++) {
 | 
|---|
 | 188 |       if (!InputFile.eof()) {
 | 
|---|
 | 189 |         InputFile >> TEFactor;
 | 
|---|
 | 190 |         (*runner).second.second = TEFactor;
 | 
|---|
 | 191 |         LOG(2, "INFO: Setting " << ++NumberOfFragments << " fragment's TEFactor to " << (*runner).second.second << ".");
 | 
|---|
 | 192 |       } else {
 | 
|---|
 | 193 |         status = false;
 | 
|---|
 | 194 |         break;
 | 
|---|
 | 195 |       }
 | 
|---|
 | 196 |     }
 | 
|---|
 | 197 |     // 4. Free and done
 | 
|---|
 | 198 |     InputFile.close();
 | 
|---|
 | 199 |     LOG(1, "INFO: done.");
 | 
|---|
 | 200 |   } else {
 | 
|---|
 | 201 |     LOG(1, "INFO: File " << filename << " not found.");
 | 
|---|
 | 202 |     status = false;
 | 
|---|
 | 203 |   }
 | 
|---|
 | 204 | 
 | 
|---|
 | 205 |   return status;
 | 
|---|
 | 206 | };
 | 
|---|
 | 207 | 
 | 
|---|
 | 208 | /** Stores TEFactors to file.
 | 
|---|
 | 209 |  * \param *out output stream for debugging
 | 
|---|
 | 210 |  * \param KeySetList Graph with factors
 | 
|---|
 | 211 |  * \param *path path to file
 | 
|---|
 | 212 |  * \return true - file written successfully, false - writing failed
 | 
|---|
 | 213 |  */
 | 
|---|
| [ca8bea] | 214 | bool Graph::StoreTEFactorsFile(char *path) const
 | 
|---|
| [75363b] | 215 | {
 | 
|---|
 | 216 |   ofstream output;
 | 
|---|
 | 217 |   bool status =  true;
 | 
|---|
 | 218 |   string line;
 | 
|---|
 | 219 | 
 | 
|---|
 | 220 |   // open TEFactors file
 | 
|---|
 | 221 |   line = path;
 | 
|---|
 | 222 |   line.append("/");
 | 
|---|
 | 223 |   line += FRAGMENTPREFIX;
 | 
|---|
 | 224 |   line += TEFACTORSFILE;
 | 
|---|
 | 225 |   output.open(line.c_str(), ios::out);
 | 
|---|
 | 226 |   LOG(1, "INFO: Saving TEFactors of the total graph ... ");
 | 
|---|
 | 227 |   if(output != NULL) {
 | 
|---|
| [ca8bea] | 228 |     for(Graph::const_iterator runner = begin(); runner != end(); runner++)
 | 
|---|
| [75363b] | 229 |       output << (*runner).second.second << endl;
 | 
|---|
 | 230 |     LOG(1, "INFO: done." << endl);
 | 
|---|
 | 231 |   } else {
 | 
|---|
 | 232 |     ELOG(2, "INFO: failed to open " << line << "." << endl);
 | 
|---|
 | 233 |     status = false;
 | 
|---|
 | 234 |   }
 | 
|---|
 | 235 |   output.close();
 | 
|---|
 | 236 | 
 | 
|---|
 | 237 |   return status;
 | 
|---|
 | 238 | };
 | 
|---|
 | 239 | 
 | 
|---|
 | 240 | /** For a given graph, sorts KeySets into a (index, keyset) map.
 | 
|---|
 | 241 |  * \return ref to allocated map from index to keyset
 | 
|---|
 | 242 |  */
 | 
|---|
| [730d7a] | 243 | AdaptivityMap * Graph::GraphToAdaptivityMap() const
 | 
|---|
| [75363b] | 244 | {
 | 
|---|
| [730d7a] | 245 |   AdaptivityMap *IndexKeySetList = new AdaptivityMap;
 | 
|---|
| [75363b] | 246 |   for(const_iterator runner = begin(); runner != end(); runner++) {
 | 
|---|
 | 247 |     IndexKeySetList->insert( pair<int,KeySet>(runner->second.first,runner->first) );
 | 
|---|
 | 248 |   }
 | 
|---|
 | 249 |   return IndexKeySetList;
 | 
|---|
 | 250 | };
 | 
|---|
| [bfbd4a] | 251 | 
 | 
|---|
 | 252 | /** For a given molecule, returns each keyset local to this molecule.
 | 
|---|
 | 253 |  *
 | 
|---|
 | 254 |  * @param mol desired molecule
 | 
|---|
 | 255 |  * @return graph with local keysets (local ids!)
 | 
|---|
 | 256 |  */
 | 
|---|
 | 257 | Graph Graph::getLocalGraph(const molecule* mol) const
 | 
|---|
 | 258 | {
 | 
|---|
 | 259 |   Graph LocalGraph;
 | 
|---|
 | 260 |   // go through all key sets, assuming ids are global
 | 
|---|
 | 261 |   for (const_iterator iter = begin(); iter != end(); ++iter) {
 | 
|---|
 | 262 |     // check whether each keyset's GLOBAL id is contained in mol_ids (via conversion through World)
 | 
|---|
 | 263 |     KeySet temp;
 | 
|---|
 | 264 |     for (KeySet::const_iterator keyiter = iter->first.begin();
 | 
|---|
 | 265 |         keyiter != iter->first.end(); ++keyiter) {
 | 
|---|
 | 266 |       const size_t globalid = *keyiter;
 | 
|---|
 | 267 |       const atom *Walker = World::getInstance().getAtom(AtomById(globalid));
 | 
|---|
 | 268 |       if (Walker != NULL) {
 | 
|---|
 | 269 |         if (Walker->getMolecule() != mol) {
 | 
|---|
 | 270 |           break;
 | 
|---|
 | 271 |         } else {
 | 
|---|
 | 272 |           // and store the LOCAL number
 | 
|---|
 | 273 |           const size_t localid = Walker->getNr();
 | 
|---|
 | 274 |           temp.insert(localid);
 | 
|---|
 | 275 |         }
 | 
|---|
 | 276 |       } else {
 | 
|---|
 | 277 |         ELOG(0, "Id " << globalid << " is not associated with any atom.");
 | 
|---|
 | 278 |         break;
 | 
|---|
 | 279 |       }
 | 
|---|
 | 280 |     }
 | 
|---|
 | 281 |     if (temp.size() == iter->first.size()) {
 | 
|---|
 | 282 |       // if so, add to LocalGraph
 | 
|---|
 | 283 |       LocalGraph.insert( std::make_pair(temp, iter->second) );
 | 
|---|
 | 284 |     }
 | 
|---|
 | 285 |   }
 | 
|---|
 | 286 | 
 | 
|---|
 | 287 |   return LocalGraph;
 | 
|---|
 | 288 | }
 | 
|---|