source: src/Actions/MoleculeAction/CreateMicelleAction.cpp@ 42127c

/*
 * Project: MoleCuilder
 * Description: creates and alters molecular systems
 * Copyright (C)  2010 University of Bonn. All rights reserved.
 * Please see the LICENSE file or "Copyright notice" in builder.cpp for details.
 */

/*
 * CreateMicelleAction.cpp
 *
 *  Created on: Sept 29, 2010
 *      Author: dueck
 */

// include config.h
#ifdef HAVE_CONFIG_H
#include <config.h>
#endif
#include "CodePatterns/MemDebug.hpp"

#include "Actions/ActionHistory.hpp"
#include "Actions/ActionRegistry.hpp"
#include "Actions/GraphAction/SubgraphDissectionAction.hpp"
#include "Actions/MoleculeAction/RotateToPrincipalAxisSystemAction.hpp"
#include "Descriptors/AtomIdDescriptor.hpp"
#include "Descriptors/MoleculeDescriptor.hpp"

#include "atom.hpp"
#include "Bond/bond.hpp"
#include "CodePatterns/Assert.hpp"
#include "CodePatterns/Log.hpp"
#include "CodePatterns/Verbose.hpp"
#include "LinearAlgebra/Line.hpp"
#include "molecule.hpp"
#include "MoleculeListClass.hpp"
#include "World.hpp"

#include <iostream>
#include <string>
#include <vector>

#include "Parser/PdbParser.hpp"
#include "Parser/TremoloParser.hpp"
#include "Parser/XyzParser.hpp"
#include "Parser/FormatParserStorage.hpp"
#include "Shapes/BaseShapes.hpp"
#include "Shapes/ShapeOps.hpp"

#include "Actions/MoleculeAction/CreateMicelleAction.hpp"

#include "CreateMicelleAction.def"
#include "Action_impl_pre.hpp"

#include "UIElements/UIFactory.hpp"
#include "UIElements/Dialog.hpp"
#include "Actions/ValueStorage.hpp"

using namespace MoleCuilder;

/** =========== define the function ====================== */
Action::state_ptr MoleculeCreateMicelleAction::performCall()
{
  typedef std::vector <atom *> AtomVector;
  typedef std::vector <molecule *> MoleculeVector;

  getParametersfromValueStorage();
        
  AtomVector ever = World::getInstance().getAllAtoms();

  // as all parsed atoms go into same molecule
  // we don't need to create one and add them all to it
  MoleculeVector all = World::getInstance().getSelectedMolecules();

  // return success when none are selected
  if (all.empty())
    return Action::success;

  // check that there is precisely one molecule
  ASSERT(all.size() == 1,
      "MoleculeCreateMicelleAction::performCall() - more than one molecule selected, this is not supported.");

  molecule *stick = *(all.begin());


  // create undo state for original molecule
  std::vector<Vector> stickPositions;
  for (molecule::const_iterator iter = stick->begin(); iter != stick->end(); ++iter)
    stickPositions.push_back((*iter)->getPosition());

  // determine principal axis system and rotate such that greatest extension is along z axis
  if(stick->size() > 1) { // only rotate if more than one molecule present
    Vector den(0.0,0.0,1.0);

    MoleculeRotateToPrincipalAxisSystem(den);
  }

  // center molecule
  stick->CenterOrigin();

  /// Align molecule with its PAS multiple times with the some surface

  // get points on surface

  //double radius= 1.5*sqrt(pow(1.55, 2)*params.N);

  Shape s = resize(Sphere(), params.radius);
  std::vector<Vector> pt = s.getHomogeneousPointsOnSurface(params.N);

  // mirror along x-y plane

  // additionally invert molecule if desired
  if (params.DoRotate) {
    for (molecule::iterator it2=stick->begin();it2 !=stick->end();++it2)
    {
      Vector pos= (**it2).getPosition();
      pos[2]= - pos[2];// -Min[2]
      (**it2).setPosition(pos);
    }
  }

  // shift molecule by its extension along z axis

  for (molecule::iterator it2=stick->begin();it2 !=stick->end();++it2)
  {
          Vector pos= (**it2).getPosition();
          pos[2]=pos[2]+params.radius;
          (**it2).setPosition(pos);
  }

  // copy molecule N times and rotate it to point radially away from surface

  //double MYEPSILON=1e-10;

  // gather created molecules for undo state
  std::vector<molecule *> CreatedMolecules;

  // TODO: remove this when we do not need MoleculeListClass anymore
  MoleculeListClass *&molecules = World::getInstance().getMolecules();

  size_t ka = 0;
  for (; ka<pt.size()-1; ka++)
  {
    LOG(1, "STATUS: Creating " << ka+1 << " copy of tenside molecule 'stick' with " << stick->getAtomCount() << " atoms, ");
    molecule *Tensid=stick->CopyMolecule();
    molecules->insert(Tensid);
    CreatedMolecules.push_back(Tensid);

    Vector ZAxis(Vector(0.0,0.0,1.0));
    Vector Axis(pt[ka]);
    const double alpha = ZAxis.Angle(Axis);
    Axis.VectorProduct(ZAxis);
    Line RotationAxis(Vector(0.0,0.0,1.0), Axis);       // pt is the current Vector of point on surface
    LOG(2, "INFO: Rotating along axis " << RotationAxis << " by angle " << alpha << ".");

    for (molecule::iterator it2=Tensid->begin();it2 !=Tensid->end();++it2)
    {
      (*it2)->setPosition(RotationAxis.rotateVector((*it2)->getPosition(),alpha));
      *(*it2)+=params.center;
    }

    Tensid=NULL;
  }

  // shift molecule at given position on surface
  {
    LOG(1, "STATUS: Shifting " << ka+1 << " copy of tenside molecule 'stick' with " << stick->getAtomCount() << " atoms, ");
    molecule *Tensid=stick;
    // this molecule was present before, hence do not remove it!
    //CreatedMolecules.push_back(Tensid);

    Vector ZAxis(Vector(0.0,0.0,1.0));
    Vector Axis(pt[pt.size()-1]);
    const double alpha = ZAxis.Angle(Axis);
    Axis.VectorProduct(ZAxis);
    Line RotationAxis(Vector(0.0,0.0,1.0), Axis);       // pt is the current Vector of point on surface
    LOG(2, "INFO: Rotating along axis " << RotationAxis << " by angle " << alpha << ".");

    for (molecule::iterator it2=Tensid->begin();it2 !=Tensid->end();++it2)
    {
      (*it2)->setPosition(RotationAxis.rotateVector((*it2)->getPosition(),alpha));
      *(*it2)+=params.center;
    }
  }


  return Action::state_ptr(new MoleculeCreateMicelleState(CreatedMolecules, stick, stickPositions, params));
}

Action::state_ptr MoleculeCreateMicelleAction::performUndo(Action::state_ptr _state)
{
  MoleculeCreateMicelleState *state = assert_cast<MoleculeCreateMicelleState*>(_state.get());

  // removing all created molecules
  MoleculeListClass *&molecules = World::getInstance().getMolecules();
  LOG(0, "STATUS: Removing " << state->CreatedMolecules.size() << " molecules again.");
  for(std::vector<molecule *>::const_iterator iter = state->CreatedMolecules.begin();
     iter != state->CreatedMolecules.end();
     ++iter) {
    LOG(2, "INFO: Removing " << (*iter)->getName());
    for(molecule::iterator AtomRunner = (*iter)->begin();
      !(*iter)->empty();
      AtomRunner = (*iter)->begin()) {
      World::getInstance().destroyAtom(*AtomRunner);
    }
    //LOG(3, "INFO: molecules contains still " << molecules->ListOfMolecules.size() << " molecules.");
    molecules->erase(*iter);
    World::getInstance().destroyMolecule(*iter);
  }
  
  // undo changes to original molecule
  ASSERT(state->stickPositions.size() == state->stick->size(),
      "MoleculeCreateMicelleAction::performUndo() - number of stored positions ("
      +toString(state->stickPositions.size())+") and number of atoms ("
      +toString(state->stick->size())+") do not match.");
  std::vector<Vector>::const_iterator PosIter = state->stickPositions.begin();
  molecule::iterator iter = state->stick->begin();
  for (; iter != state->stick->end(); ++iter, ++PosIter)
    (*iter)->setPosition(*PosIter);

  return Action::state_ptr(_state);
}

Action::state_ptr MoleculeCreateMicelleAction::performRedo(Action::state_ptr _state)
{
  LOG(0, "STATUS: CreateMicelleAction::performRedo() is nit implemented yet.");
  return Action::failure;
}


bool MoleculeCreateMicelleAction::canUndo()
{
  return true;
}

bool MoleculeCreateMicelleAction::shouldUndo()
{
  return true;
}

/** =========== end of function ====================== */