source: src/Fragmentation/Exporters/ExportGraph_ToFiles.cpp@ cee9e8

/*
 * Project: MoleCuilder
 * Description: creates and alters molecular systems
 * Copyright (C)  2012 University of Bonn. All rights reserved.
 * Copyright (C)  2013 Frederik Heber. 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 <http://www.gnu.org/licenses/>.
 */

/*
 * ExportGraph_ToFiles.cpp
 *
 *  Created on: 08.03.2012
 *      Author: heber
 */

// include config.h
#ifdef HAVE_CONFIG_H
#include <config.h>
#endif

#include "CodePatterns/MemDebug.hpp"

#include "ExportGraph_ToFiles.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 "Parser/FormatParserStorage.hpp"
#include "World.hpp"

/** Constructor for class ExportGraph_ToFiles.
 *
 * @param _graph instance of Graph containing keyset of each fragment
 * @param _treatment whether to always add already present hydrogens or not
 * @param _saturation whether to saturate dangling bonds with hydrogen or not
 * @param _globalsaturationpositions possibly empty map with global information
 *        where to place saturation hydrogens to fulfill consistency principle
 */
ExportGraph_ToFiles::ExportGraph_ToFiles(
    const Graph &_graph,
    const enum HydrogenTreatment _treatment,
    const enum HydrogenSaturation _saturation,
    const SaturatedFragment::GlobalSaturationPositions_t &_globalsaturationpositions) :
                ExportGraph(_graph, _treatment, _saturation, _globalsaturationpositions)
{}

/** Destructor of class ExportGraph_ToFiles.
 *
 * We free all created molecules again and also removed their copied atoms.
 */
ExportGraph_ToFiles::~ExportGraph_ToFiles()
{}

/** Returns a string with \a i prefixed with 0s to match order of total number of molecules in digits.
 * \param FragmentNumber total number of fragments to determine necessary number of digits
 * \param digits number to create with 0 prefixed
 * \return allocated(!) char array with number in digits, ten base.
 */
static char *FixedDigitNumber(const int FragmentNumber, const int digits)
{
  char *returnstring;
  int number = FragmentNumber;
  int order = 0;
  while (number != 0) { // determine number of digits needed
    number = (int)floor(((double)number / 10.));
    order++;
    //LOG(0, "Number is " << number << ", order is " << order << ".");
  }
  // allocate string
  returnstring = new char[order + 2];
  // terminate  and fill string array from end backward
  returnstring[order] = '\0';
  number = digits;
  for (int i=order;i--;) {
    returnstring[i] = '0' + (char)(number % 10);
    number = (int)floor(((double)number / 10.));
  }
  //LOG(0, returnstring);
  return returnstring;
};

/** Actual implementation of the export to files function.
 */
void ExportGraph_ToFiles::operator()()
{
  LOG(1, "INFO: Writing " << TotalGraph.size() << " possible bond fragmentation configs.");
  size_t FragmentCounter = 0;
  char *FragmentNumber = NULL;
  string filename(prefix);
  filename += FORCESFILE;
  std::ofstream ForcesFile(filename.c_str());
  SortIndex_t SortIndex;

  // ===== 9. Save fragments' configuration and keyset files et al to disk ===
  bool write_status = true;
  ExportGraph::SaturatedFragment_ptr CurrentFragment = getNextFragment();
  for (; (CurrentFragment != NULL) && (CurrentFragment->getKeySet() != ExportGraph::EmptySet);
      CurrentFragment = getNextFragment()) {
    const KeySet &set = CurrentFragment->getKeySet();
    LOG(2, "INFO: Writing bond fragments for set " << set << ".");
    // store config in stream
    {
      // open file
      FragmentNumber = FixedDigitNumber(TotalGraph.size(), FragmentCounter++);
      storeFragmentForAllTypes(
          CurrentFragment, FragmentNumber, FragmentCounter-1);
      delete[](FragmentNumber);
    }
    // store force index reference file
    write_status = write_status
        && appendToForcesFile(CurrentFragment, ForcesFile, SortIndex);
    // explicitly release fragment
    CurrentFragment.reset();
  }
  if (CurrentFragment == NULL) {
    ELOG(1, "Some error while obtaining the next fragment occured.");
    return;
  }
  ForcesFile.close();

  if (write_status)
    LOG(1, "All configs written.");
  else
    LOG(1, "Some config writing failed.");

  // store keysets file
  TotalGraph.StoreKeySetFile(prefix);

  // store Hydrogen saturation correction file
//  BondFragments.AddHydrogenCorrection(prefix);

  // restore orbital and Stop values
  //CalculateOrbitals(*configuration);
}

bool ExportGraph_ToFiles::storeFragmentForAllTypes(
    SaturatedFragment_ptr &CurrentFragment,
    char *FragmentNumber,
    size_t FragmentCounter) const
{
  bool write_status = true;

  // go through all desired types
  for (std::vector<std::string>::const_iterator typeiter = typelist.begin();
      typeiter != typelist.end(); ++typeiter) {
    const std::string &typeName = *typeiter;
    const ParserTypes type =
        FormatParserStorage::getInstance().getTypeFromName(typeName);
    // create filenname and open
    const std::string FragmentName =
        prefix + FragmentNumber + "." + FormatParserStorage::getInstance().getSuffixFromType(type);
    std::ofstream outputFragment(FragmentName.c_str(), ios::out);

    // write to this stream
    {
      std::stringstream output;
      output << "INFO: Saving bond fragment No. " << FragmentNumber << "/"
          << FragmentCounter << " as " << typeName << " ... ";
      const bool intermediate_result = CurrentFragment->OutputConfig(outputFragment,type);
      write_status &= intermediate_result;
      if (intermediate_result)
        output << " done.";
      else
        output << " failed.";
      LOG(2, output.str());
    }

    // close file
    outputFragment.close();
    outputFragment.clear();
  }

  return write_status;
}

bool ExportGraph_ToFiles::appendToForcesFile(
    SaturatedFragment_ptr &CurrentFragment,
    std::ostream &ForcesFile,
    const SortIndex_t &SortIndex) const
{
  bool status = true;
//  periodentafel *periode=World::getInstance().getPeriode();

  // open file for the force factors
  if (ForcesFile.good()) {
    //output << prefix << "Forces" << endl;
    const KeySet &FullMolecule = CurrentFragment->getFullMolecule();
    const KeySet &SaturationHydrogens = CurrentFragment->getSaturationHydrogens();
    for (KeySet::const_iterator keyiter = FullMolecule.begin();
        keyiter != FullMolecule.end();
        ++keyiter) {
      if (SaturationHydrogens.find(*keyiter) == SaturationHydrogens.end()) {
        ForcesFile << SortIndex.find(*keyiter) << "\t";
      } else {
        // otherwise a -1 to indicate an added saturation hydrogen
        ForcesFile << "-1\t";
      }
    }
//    for (MoleculeList::iterator ListRunner = ListOfMolecules.begin(); ListRunner != ListOfMolecules.end(); ListRunner++) {
//      periodentafel::const_iterator elemIter;
//      for(elemIter=periode->begin();elemIter!=periode->end();++elemIter){
//        if ((*ListRunner)->hasElement((*elemIter).first)) { // if this element got atoms
//          for(molecule::iterator atomIter = (*ListRunner)->begin(); atomIter !=(*ListRunner)->end();++atomIter){
//            if ((*atomIter)->getType()->getAtomicNumber() == (*elemIter).first) {
//              if (((*atomIter)->GetTrueFather() != NULL) && ((*atomIter)->GetTrueFather() != (*atomIter))) {// if there is a rea
//                const atomId_t fatherid = (*atomIter)->GetTrueFather()->getId();
//                ForcesFile << SortIndex.find(fatherid) << "\t";
//              } else
//                // otherwise a -1 to indicate an added saturation hydrogen
//                ForcesFile << "-1\t";
//            }
//          }
//        }
//      }
//      ForcesFile << endl;
//    }
    ForcesFile << std::endl;
  } else {
    status = false;
    ELOG(1, "Failure on appending to ForcesFile.");
  }

  return status;
}