Index: src/Actions/AnalysisAction/DipoleAngularCorrelationAction.cpp
===================================================================
--- src/Actions/AnalysisAction/DipoleAngularCorrelationAction.cpp	(revision a58c16f47065e9446b7916a407271c6544180022)
+++ src/Actions/AnalysisAction/DipoleAngularCorrelationAction.cpp	(revision 99db9baca1782e208fb17d2cc6ace3f4cadb568b)
@@ -64,6 +64,8 @@
 
   // get selected atoms
-  std::vector<atom*> old_atom_selection = World::getInstance().getSelectedAtoms();
-  std::vector<molecule*> old_molecule_selection = World::getInstance().getSelectedMolecules();
+  std::vector<const atom*> old_atom_selection =
+      const_cast<const World &>(World::getInstance()).getSelectedAtoms();
+  std::vector<const molecule*> old_molecule_selection =
+      const_cast<const World &>(World::getInstance()).getSelectedMolecules();
 
   // get current time step
@@ -75,5 +77,6 @@
   World::getInstance().clearMoleculeSelection(); // TODO: This should be done in setTime or where molecules are re-done
   World::getInstance().selectAllMolecules(MoleculeByFormula(DipoleFormula));
-  std::vector<molecule *> molecules = World::getInstance().getSelectedMolecules();
+  std::vector<const molecule *> molecules =
+      const_cast<const World &>(World::getInstance()).getSelectedMolecules();
   std::map<atomId_t, Vector> ZeroVector = CalculateZeroAngularDipole(molecules);
 
@@ -81,5 +84,6 @@
   World::getInstance().clearAtomSelection();
   World::getInstance().selectAllAtoms(AtomsByMoleculeSelection());
-  std::vector<atom *> atoms = World::getInstance().getSelectedAtoms();
+  std::vector<const atom *> atoms = const_cast<const World &>(World::getInstance()).
+      getSelectedAtoms();
   if (atoms.empty()) {
     STATUS("Formula "+toString(DipoleFormula)+" selects no atoms.");
@@ -130,9 +134,9 @@
   // reset to old selections
   World::getInstance().clearAtomSelection();
-  BOOST_FOREACH(atom *_atom, old_atom_selection) {
+  BOOST_FOREACH(const atom *_atom, old_atom_selection) {
     World::getInstance().selectAtom(_atom);
   }
   World::getInstance().clearMoleculeSelection();
-  BOOST_FOREACH(molecule *_mol, old_molecule_selection) {
+  BOOST_FOREACH(const molecule *_mol, old_molecule_selection) {
     World::getInstance().selectMolecule(_mol);
   }
Index: src/Actions/AnalysisAction/DipoleCorrelationAction.cpp
===================================================================
--- src/Actions/AnalysisAction/DipoleCorrelationAction.cpp	(revision a58c16f47065e9446b7916a407271c6544180022)
+++ src/Actions/AnalysisAction/DipoleCorrelationAction.cpp	(revision 99db9baca1782e208fb17d2cc6ace3f4cadb568b)
@@ -66,5 +66,6 @@
   binoutput.open(params.binoutputname.get().string().c_str());
   DipoleCorrelationMap *correlationmap = NULL;
-  std::vector<molecule*> molecules = World::getInstance().getSelectedMolecules();
+  const std::vector<const molecule*> molecules =
+      const_cast<const World &>(World::getInstance()).getSelectedMolecules();
   STATUS("There are "+toString(molecules.size())+" selected molecules.");
   ASSERT(!params.periodic.get(), "AnalysisDipoleCorrelationAction() - periodic case not implemented.");
Index: src/Actions/AnalysisAction/PointCorrelationAction.cpp
===================================================================
--- src/Actions/AnalysisAction/PointCorrelationAction.cpp	(revision a58c16f47065e9446b7916a407271c6544180022)
+++ src/Actions/AnalysisAction/PointCorrelationAction.cpp	(revision 99db9baca1782e208fb17d2cc6ace3f4cadb568b)
@@ -71,5 +71,6 @@
   for(std::vector<const element *>::const_iterator iter = params.elements.get().begin(); iter != params.elements.get().end(); ++iter)
     cout << "element is " << (*iter)->getSymbol() << endl;
-  std::vector<molecule*> molecules = World::getInstance().getSelectedMolecules();
+  const std::vector<const molecule*> molecules =
+      const_cast<const World &>(World::getInstance()).getSelectedMolecules();
   if (params.periodic.get())
     correlationmap  = PeriodicCorrelationToPoint(molecules, params.elements.get(), &params.Point.get(), ranges);
Index: src/Actions/AnalysisAction/SurfaceCorrelationAction.cpp
===================================================================
--- src/Actions/AnalysisAction/SurfaceCorrelationAction.cpp	(revision a58c16f47065e9446b7916a407271c6544180022)
+++ src/Actions/AnalysisAction/SurfaceCorrelationAction.cpp	(revision 99db9baca1782e208fb17d2cc6ace3f4cadb568b)
@@ -96,5 +96,6 @@
 
   // correlate
-  std::vector<molecule*> molecules = World::getInstance().getSelectedMolecules();
+  const std::vector<const molecule*> molecules =
+      const_cast<const World &>(World::getInstance()).getSelectedMolecules();
   std::cout << "There are " << molecules.size() << " selected molecules." << std::endl;
   CorrelationToSurfaceMap *surfacemap = NULL;
Index: src/Actions/AtomAction/RotateAroundOriginByAngleAction.cpp
===================================================================
--- src/Actions/AtomAction/RotateAroundOriginByAngleAction.cpp	(revision a58c16f47065e9446b7916a407271c6544180022)
+++ src/Actions/AtomAction/RotateAroundOriginByAngleAction.cpp	(revision 99db9baca1782e208fb17d2cc6ace3f4cadb568b)
@@ -76,5 +76,5 @@
   }
   LOG(0, "done.");
-  return ActionState::ptr(new AtomRotateAroundOriginByAngleState(World::getInstance().getSelectedAtoms(), params));
+  return ActionState::ptr(new AtomRotateAroundOriginByAngleState(selectedAtoms, params));
 }
 
Index: src/Actions/MoleculeAction/SaveTemperatureAction.cpp
===================================================================
--- src/Actions/MoleculeAction/SaveTemperatureAction.cpp	(revision a58c16f47065e9446b7916a407271c6544180022)
+++ src/Actions/MoleculeAction/SaveTemperatureAction.cpp	(revision 99db9baca1782e208fb17d2cc6ace3f4cadb568b)
@@ -60,7 +60,8 @@
   ofstream output;
   output.open(params.temperaturefile.get().string().c_str(), ios::trunc);
-  AtomSetMixin<std::vector<atom *> > set(World::getInstance().getSelectedAtoms());
+  AtomSetMixin<std::vector<const atom *> > set =
+      const_cast<const World &>(World::getInstance()).getSelectedAtoms();
   const size_t MDSteps = set.getMaxTrajectorySize();
-  OutputTemperature<std::vector<atom *> > writer(set);
+  OutputTemperature<std::vector<const atom *> > writer(set);
   if (output.fail() || !writer((ofstream * const) &output, 0, MDSteps)) {
     STATUS("File could not be written.");
Index: src/Actions/ParserAction/SaveSelectedAtomsAsExtTypesAction.cpp
===================================================================
--- src/Actions/ParserAction/SaveSelectedAtomsAsExtTypesAction.cpp	(revision a58c16f47065e9446b7916a407271c6544180022)
+++ src/Actions/ParserAction/SaveSelectedAtomsAsExtTypesAction.cpp	(revision 99db9baca1782e208fb17d2cc6ace3f4cadb568b)
@@ -65,5 +65,8 @@
     boost::filesystem::ofstream test;
     test.open(params.filename.get());
-    const bool status = parser.saveAtomsInExttypes(test, World::getInstance().getSelectedAtoms(), params.id.get());
+    const bool status = parser.saveAtomsInExttypes(
+        test,
+        const_cast<const World &>(World::getInstance()).getSelectedAtoms(),
+        params.id.get());
     test.close();
 
Index: src/Actions/SelectionAction/Atoms/AllAtomsAction.cpp
===================================================================
--- src/Actions/SelectionAction/Atoms/AllAtomsAction.cpp	(revision a58c16f47065e9446b7916a407271c6544180022)
+++ src/Actions/SelectionAction/Atoms/AllAtomsAction.cpp	(revision 99db9baca1782e208fb17d2cc6ace3f4cadb568b)
@@ -54,5 +54,6 @@
 /** =========== define the function ====================== */
 ActionState::ptr SelectionAllAtomsAction::performCall() {
-  std::vector<atom *> selectedAtoms = World::getInstance().getSelectedAtoms();
+  const std::vector<const atom *> selectedAtoms =
+      const_cast<const World &>(World::getInstance()).getSelectedAtoms();
   LOG(1, "Selecting all atoms.");
   World::getInstance().selectAllAtoms(AllAtoms());
@@ -66,5 +67,5 @@
   World::getInstance().clearAtomSelection();
 
-  BOOST_FOREACH(atom *_atom, state->selectedAtoms)
+  BOOST_FOREACH(const atom *_atom, state->selectedAtoms)
     World::getInstance().selectAtom(_atom);
 
Index: src/Actions/SelectionAction/Atoms/AllAtomsAction.def
===================================================================
--- src/Actions/SelectionAction/Atoms/AllAtomsAction.def	(revision a58c16f47065e9446b7916a407271c6544180022)
+++ src/Actions/SelectionAction/Atoms/AllAtomsAction.def	(revision 99db9baca1782e208fb17d2cc6ace3f4cadb568b)
@@ -21,5 +21,5 @@
 #undef paramvalids
 
-#define statetypes (std::vector<atom*>)
+#define statetypes (std::vector<const atom*>)
 #define statereferences (selectedAtoms)
 
Index: src/Actions/SelectionAction/Atoms/ClearAllAtomsAction.cpp
===================================================================
--- src/Actions/SelectionAction/Atoms/ClearAllAtomsAction.cpp	(revision a58c16f47065e9446b7916a407271c6544180022)
+++ src/Actions/SelectionAction/Atoms/ClearAllAtomsAction.cpp	(revision 99db9baca1782e208fb17d2cc6ace3f4cadb568b)
@@ -54,5 +54,6 @@
 /** =========== define the function ====================== */
 ActionState::ptr SelectionClearAllAtomsAction::performCall() {
-  std::vector<atom *> selectedAtoms = World::getInstance().getSelectedAtoms();
+  const std::vector<const atom *> selectedAtoms =
+      const_cast<const World &>(World::getInstance()).getSelectedAtoms();
   LOG(1, "Clearing atoms selection.");
   World::getInstance().clearAtomSelection();
@@ -64,5 +65,5 @@
 
   World::getInstance().clearAtomSelection();
-  BOOST_FOREACH(atom *_atom, state->selectedAtoms)
+  BOOST_FOREACH(const atom *_atom, state->selectedAtoms)
     World::getInstance().selectAtom(_atom);
 
Index: src/Actions/SelectionAction/Atoms/ClearAllAtomsAction.def
===================================================================
--- src/Actions/SelectionAction/Atoms/ClearAllAtomsAction.def	(revision a58c16f47065e9446b7916a407271c6544180022)
+++ src/Actions/SelectionAction/Atoms/ClearAllAtomsAction.def	(revision 99db9baca1782e208fb17d2cc6ace3f4cadb568b)
@@ -20,5 +20,5 @@
 #undef paramvalids
 
-#define statetypes (std::vector<atom*>)
+#define statetypes (std::vector<const atom*>)
 #define statereferences (selectedAtoms)
 
Index: src/Actions/SelectionAction/Atoms/NotAllAtomsAction.cpp
===================================================================
--- src/Actions/SelectionAction/Atoms/NotAllAtomsAction.cpp	(revision a58c16f47065e9446b7916a407271c6544180022)
+++ src/Actions/SelectionAction/Atoms/NotAllAtomsAction.cpp	(revision 99db9baca1782e208fb17d2cc6ace3f4cadb568b)
@@ -54,5 +54,6 @@
 /** =========== define the function ====================== */
 ActionState::ptr SelectionNotAllAtomsAction::performCall() {
-  std::vector<atom *> selectedAtoms = World::getInstance().getSelectedAtoms();
+  const std::vector<const atom *> selectedAtoms =
+      const_cast<const World &>(World::getInstance()).getSelectedAtoms();
   LOG(1, "Unselecting all atoms.");
   World::getInstance().clearAtomSelection();
@@ -65,5 +66,5 @@
 
   World::getInstance().clearAtomSelection();
-  BOOST_FOREACH(atom *_atom, state->unselectedAtoms)
+  BOOST_FOREACH(const atom *_atom, state->unselectedAtoms)
     World::getInstance().selectAtom(_atom);
 
Index: src/Actions/SelectionAction/Atoms/NotAllAtomsAction.def
===================================================================
--- src/Actions/SelectionAction/Atoms/NotAllAtomsAction.def	(revision a58c16f47065e9446b7916a407271c6544180022)
+++ src/Actions/SelectionAction/Atoms/NotAllAtomsAction.def	(revision 99db9baca1782e208fb17d2cc6ace3f4cadb568b)
@@ -21,5 +21,5 @@
 #undef paramvalids
 
-#define statetypes (std::vector<atom*>)
+#define statetypes (std::vector<const atom*>)
 #define statereferences (unselectedAtoms)
 
Index: src/Actions/SelectionAction/Atoms/PopAtomsAction.cpp
===================================================================
--- src/Actions/SelectionAction/Atoms/PopAtomsAction.cpp	(revision a58c16f47065e9446b7916a407271c6544180022)
+++ src/Actions/SelectionAction/Atoms/PopAtomsAction.cpp	(revision 99db9baca1782e208fb17d2cc6ace3f4cadb568b)
@@ -50,13 +50,9 @@
 ActionState::ptr SelectionPopAtomsAction::performCall() {
   // create undo state
-  const std::vector<atom *> selected_atoms = World::getInstance().getSelectedAtoms();
-  std::vector<atomId_t> selected_atomids(selected_atoms.size(), (atomId_t)-1);
-  std::transform(
-      selected_atoms.begin(), selected_atoms.end(),
-      selected_atomids.begin(),
-      boost::bind(&atom::getId, _1));
+  const std::vector<const atom *> selectedAtoms =
+      const_cast<const World &>(World::getInstance()).getSelectedAtoms();
   LOG(1, "Popping atom selection.");
   World::getInstance().popAtomSelection();
-  return ActionState::ptr(new SelectionPopAtomsState(selected_atomids, params));
+  return ActionState::ptr(new SelectionPopAtomsState(selectedAtoms, params));
 }
 
@@ -67,7 +63,7 @@
 
   // and re-create present selection
-  for(std::vector<atomId_t>::const_iterator iter = state->selected_atomids.begin();
-      iter != state->selected_atomids.end(); ++iter)
-    World::getInstance().selectAtom(*iter);
+  World::getInstance().clearAtomSelection();
+  BOOST_FOREACH(const atom *_atom, state->selectedAtoms)
+    World::getInstance().selectAtom(_atom);
 
   return ActionState::ptr(_state);
Index: src/Actions/SelectionAction/Atoms/PopAtomsAction.def
===================================================================
--- src/Actions/SelectionAction/Atoms/PopAtomsAction.def	(revision a58c16f47065e9446b7916a407271c6544180022)
+++ src/Actions/SelectionAction/Atoms/PopAtomsAction.def	(revision 99db9baca1782e208fb17d2cc6ace3f4cadb568b)
@@ -20,6 +20,6 @@
 #undef paramvalids
 
-#define statetypes (std::vector<atomId_t>)
-#define statereferences (selected_atomids)
+#define statetypes (std::vector<const atom *>)
+#define statereferences (selectedAtoms)
 
 // some defines for all the names, you may use ACTION, STATE and PARAMS
Index: src/Actions/SelectionAction/Molecules/AllMoleculesAction.cpp
===================================================================
--- src/Actions/SelectionAction/Molecules/AllMoleculesAction.cpp	(revision a58c16f47065e9446b7916a407271c6544180022)
+++ src/Actions/SelectionAction/Molecules/AllMoleculesAction.cpp	(revision 99db9baca1782e208fb17d2cc6ace3f4cadb568b)
@@ -55,5 +55,6 @@
 /** =========== define the function ====================== */
 ActionState::ptr SelectionAllMoleculesAction::performCall() {
-  std::vector<molecule *> selectedMolecules = World::getInstance().getSelectedMolecules();
+  const std::vector<const molecule *> selectedMolecules =
+      const_cast<const World &>(World::getInstance()).getSelectedMolecules();
   LOG(1, "Selecting all molecules.");
   World::getInstance().selectAllMolecules(AllMolecules());
@@ -66,5 +67,5 @@
 
   World::getInstance().clearMoleculeSelection();
-  BOOST_FOREACH(molecule *_mol, state->selectedMolecules)
+  BOOST_FOREACH(const molecule *_mol, state->selectedMolecules)
     World::getInstance().selectMolecule(_mol);
 
Index: src/Actions/SelectionAction/Molecules/AllMoleculesAction.def
===================================================================
--- src/Actions/SelectionAction/Molecules/AllMoleculesAction.def	(revision a58c16f47065e9446b7916a407271c6544180022)
+++ src/Actions/SelectionAction/Molecules/AllMoleculesAction.def	(revision 99db9baca1782e208fb17d2cc6ace3f4cadb568b)
@@ -20,5 +20,5 @@
 #undef paramvalids
 
-#define statetypes (std::vector<molecule*>)
+#define statetypes (std::vector<const molecule*>)
 #define statereferences (selectedMolecules)
 
Index: src/Actions/SelectionAction/Molecules/ClearAllMoleculesAction.cpp
===================================================================
--- src/Actions/SelectionAction/Molecules/ClearAllMoleculesAction.cpp	(revision a58c16f47065e9446b7916a407271c6544180022)
+++ src/Actions/SelectionAction/Molecules/ClearAllMoleculesAction.cpp	(revision 99db9baca1782e208fb17d2cc6ace3f4cadb568b)
@@ -54,5 +54,6 @@
 /** =========== define the function ====================== */
 ActionState::ptr SelectionClearAllMoleculesAction::performCall() {
-  std::vector<molecule *> selectedMolecules = World::getInstance().getSelectedMolecules();
+  const std::vector<const molecule *> selectedMolecules =
+      const_cast<const World &>(World::getInstance()).getSelectedMolecules();
   LOG(1, "Clearing molecule selection.");
   World::getInstance().clearMoleculeSelection();
@@ -65,5 +66,5 @@
 
   World::getInstance().clearMoleculeSelection();
-  BOOST_FOREACH(molecule *_mol, state->selectedMolecules)
+  BOOST_FOREACH(const molecule *_mol, state->selectedMolecules)
     World::getInstance().selectMolecule(_mol);
 
Index: src/Actions/SelectionAction/Molecules/ClearAllMoleculesAction.def
===================================================================
--- src/Actions/SelectionAction/Molecules/ClearAllMoleculesAction.def	(revision a58c16f47065e9446b7916a407271c6544180022)
+++ src/Actions/SelectionAction/Molecules/ClearAllMoleculesAction.def	(revision 99db9baca1782e208fb17d2cc6ace3f4cadb568b)
@@ -20,5 +20,5 @@
 #undef paramvalids
 
-#define statetypes (std::vector<molecule*>)
+#define statetypes (std::vector<const molecule*>)
 #define statereferences (selectedMolecules)
 
Index: src/Actions/SelectionAction/Molecules/MoleculeOfAtomAction.cpp
===================================================================
--- src/Actions/SelectionAction/Molecules/MoleculeOfAtomAction.cpp	(revision a58c16f47065e9446b7916a407271c6544180022)
+++ src/Actions/SelectionAction/Molecules/MoleculeOfAtomAction.cpp	(revision 99db9baca1782e208fb17d2cc6ace3f4cadb568b)
@@ -54,5 +54,6 @@
 /** =========== define the function ====================== */
 ActionState::ptr SelectionMoleculeOfAtomAction::performCall() {
-  std::vector<molecule *> selectedMolecules = World::getInstance().getSelectedMolecules();
+  const std::vector<const molecule *> selectedMolecules =
+      const_cast<const World &>(World::getInstance()).getSelectedMolecules();
   LOG(1, "Selecting molecule to currently selected atoms: " << World::getInstance().countSelectedAtoms());
   World::getInstance().selectAllMolecules(MoleculesByAtomSelection());
@@ -65,5 +66,5 @@
 
   World::getInstance().clearMoleculeSelection();
-  BOOST_FOREACH(molecule *_mol, state->selectedMolecules)
+  BOOST_FOREACH(const molecule *_mol, state->selectedMolecules)
     World::getInstance().selectMolecule(_mol);
 
Index: src/Actions/SelectionAction/Molecules/MoleculeOfAtomAction.def
===================================================================
--- src/Actions/SelectionAction/Molecules/MoleculeOfAtomAction.def	(revision a58c16f47065e9446b7916a407271c6544180022)
+++ src/Actions/SelectionAction/Molecules/MoleculeOfAtomAction.def	(revision 99db9baca1782e208fb17d2cc6ace3f4cadb568b)
@@ -20,5 +20,5 @@
 #undef paramvalids
 
-#define statetypes (std::vector<molecule*>)
+#define statetypes (std::vector<const molecule*>)
 #define statereferences (selectedMolecules)
 
Index: src/Actions/SelectionAction/Molecules/NotAllMoleculesAction.cpp
===================================================================
--- src/Actions/SelectionAction/Molecules/NotAllMoleculesAction.cpp	(revision a58c16f47065e9446b7916a407271c6544180022)
+++ src/Actions/SelectionAction/Molecules/NotAllMoleculesAction.cpp	(revision 99db9baca1782e208fb17d2cc6ace3f4cadb568b)
@@ -54,5 +54,6 @@
 /** =========== define the function ====================== */
 ActionState::ptr SelectionNotAllMoleculesAction::performCall() {
-  std::vector<molecule *> selectedMolecules = World::getInstance().getSelectedMolecules();
+  const std::vector<const molecule *> selectedMolecules =
+      const_cast<const World &>(World::getInstance()).getSelectedMolecules();
   LOG(1, "Unselecting all molecules.");
   World::getInstance().clearMoleculeSelection();
@@ -65,5 +66,5 @@
 
   World::getInstance().clearMoleculeSelection();
-  BOOST_FOREACH(molecule *_mol, state->unselectedMolecules)
+  BOOST_FOREACH(const molecule *_mol, state->selectedMolecules)
     World::getInstance().selectMolecule(_mol);
 
Index: src/Actions/SelectionAction/Molecules/NotAllMoleculesAction.def
===================================================================
--- src/Actions/SelectionAction/Molecules/NotAllMoleculesAction.def	(revision a58c16f47065e9446b7916a407271c6544180022)
+++ src/Actions/SelectionAction/Molecules/NotAllMoleculesAction.def	(revision 99db9baca1782e208fb17d2cc6ace3f4cadb568b)
@@ -20,6 +20,6 @@
 #undef paramvalids
 
-#define statetypes (std::vector<molecule*>)
-#define statereferences (unselectedMolecules)
+#define statetypes (std::vector<const molecule*>)
+#define statereferences (selectedMolecules)
 
 // some defines for all the names, you may use ACTION, STATE and PARAMS
Index: src/Actions/SelectionAction/Molecules/NotMoleculeOfAtomAction.cpp
===================================================================
--- src/Actions/SelectionAction/Molecules/NotMoleculeOfAtomAction.cpp	(revision a58c16f47065e9446b7916a407271c6544180022)
+++ src/Actions/SelectionAction/Molecules/NotMoleculeOfAtomAction.cpp	(revision 99db9baca1782e208fb17d2cc6ace3f4cadb568b)
@@ -54,9 +54,10 @@
 /** =========== define the function ====================== */
 ActionState::ptr SelectionNotMoleculeOfAtomAction::performCall() {
-  std::vector<molecule *> unselectedMolecules = World::getInstance().getSelectedMolecules();
+  const std::vector<const molecule *> selectedMolecules =
+      const_cast<const World &>(World::getInstance()).getSelectedMolecules();
   LOG(1, "Unselecting molecule to currently selected atoms.");
   World::getInstance().unselectAllMolecules(MoleculesByAtomSelection());
   LOG(0, World::getInstance().countSelectedMolecules() << " molecules remain selected.");
-  return ActionState::ptr(new SelectionNotMoleculeOfAtomState(unselectedMolecules, params));
+  return ActionState::ptr(new SelectionNotMoleculeOfAtomState(selectedMolecules, params));
 }
 
@@ -65,5 +66,5 @@
 
   World::getInstance().clearMoleculeSelection();
-  BOOST_FOREACH(molecule *_mol, state->unselectedMolecules)
+  BOOST_FOREACH(const molecule *_mol, state->selectedMolecules)
     World::getInstance().selectMolecule(_mol);
 
Index: src/Actions/SelectionAction/Molecules/NotMoleculeOfAtomAction.def
===================================================================
--- src/Actions/SelectionAction/Molecules/NotMoleculeOfAtomAction.def	(revision a58c16f47065e9446b7916a407271c6544180022)
+++ src/Actions/SelectionAction/Molecules/NotMoleculeOfAtomAction.def	(revision 99db9baca1782e208fb17d2cc6ace3f4cadb568b)
@@ -20,6 +20,6 @@
 #undef paramvalids
 
-#define statetypes (std::vector<molecule*>)
-#define statereferences (unselectedMolecules)
+#define statetypes (std::vector<const molecule*>)
+#define statereferences (selectedMolecules)
 
 // some defines for all the names, you may use ACTION, STATE and PARAMS
Index: src/Actions/SelectionAction/Molecules/PopMoleculesAction.cpp
===================================================================
--- src/Actions/SelectionAction/Molecules/PopMoleculesAction.cpp	(revision a58c16f47065e9446b7916a407271c6544180022)
+++ src/Actions/SelectionAction/Molecules/PopMoleculesAction.cpp	(revision 99db9baca1782e208fb17d2cc6ace3f4cadb568b)
@@ -50,13 +50,9 @@
 ActionState::ptr SelectionPopMoleculesAction::performCall() {
   // create undo state
-  const std::vector<molecule *> selected_mols = World::getInstance().getSelectedMolecules();
-  std::vector<moleculeId_t> selected_molids(selected_mols.size(), (moleculeId_t)-1);
-  std::transform(
-      selected_mols.begin(), selected_mols.end(),
-      selected_molids.begin(),
-      boost::bind(&molecule::getId, _1));
+  const std::vector<const molecule *> selectedMolecules =
+      const_cast<const World &>(World::getInstance()).getSelectedMolecules();
   LOG(1, "Popping molecule selection.");
   World::getInstance().popMoleculeSelection();
-  return ActionState::ptr(new SelectionPopMoleculesState(selected_molids, params));
+  return ActionState::ptr(new SelectionPopMoleculesState(selectedMolecules, params));
 }
 
@@ -66,7 +62,7 @@
   World::getInstance().pushMoleculeSelection();
   // re-create old selection
-  for(std::vector<moleculeId_t>::const_iterator iter = state->selected_molids.begin();
-      iter != state->selected_molids.end(); ++iter)
-    World::getInstance().selectMolecule(*iter);
+  World::getInstance().clearMoleculeSelection();
+  BOOST_FOREACH(const molecule *_mol, state->selectedMolecules)
+    World::getInstance().selectMolecule(_mol);
 
   return ActionState::ptr(_state);
Index: src/Actions/SelectionAction/Molecules/PopMoleculesAction.def
===================================================================
--- src/Actions/SelectionAction/Molecules/PopMoleculesAction.def	(revision a58c16f47065e9446b7916a407271c6544180022)
+++ src/Actions/SelectionAction/Molecules/PopMoleculesAction.def	(revision 99db9baca1782e208fb17d2cc6ace3f4cadb568b)
@@ -20,6 +20,6 @@
 #undef paramvalids
 
-#define statetypes (std::vector<moleculeId_t>)
-#define statereferences (selected_molids)
+#define statetypes (std::vector<const molecule *>)
+#define statereferences (selectedMolecules)
 
 // some defines for all the names, you may use ACTION, STATE and PARAMS
Index: src/Analysis/analysis_correlation.cpp
===================================================================
--- src/Analysis/analysis_correlation.cpp	(revision a58c16f47065e9446b7916a407271c6544180022)
+++ src/Analysis/analysis_correlation.cpp	(revision 99db9baca1782e208fb17d2cc6ace3f4cadb568b)
@@ -125,5 +125,5 @@
  * \return range with [min, max]
  */
-range<size_t> getMaximumTrajectoryBounds(const std::vector<atom *> &atoms)
+range<size_t> getMaximumTrajectoryBounds(const std::vector<const atom *> &atoms)
 {
   // get highest trajectory size
@@ -133,5 +133,5 @@
   size_t max_timesteps = std::numeric_limits<size_t>::min();
   size_t min_timesteps = std::numeric_limits<size_t>::max();
-  BOOST_FOREACH(atom *_atom, atoms) {
+  BOOST_FOREACH(const atom *_atom, atoms) {
     if (_atom->getTrajectorySize() > max_timesteps)
       max_timesteps  = _atom->getTrajectorySize();
@@ -149,17 +149,12 @@
  * \return map with orientation vector for each atomic id given in \a atoms.
  */
-std::map<atomId_t, Vector> CalculateZeroAngularDipole(const std::vector<molecule *> &molecules)
+std::map<atomId_t, Vector> CalculateZeroAngularDipole(const std::vector<const molecule *> &molecules)
 {
   // get zero orientation for each molecule.
   LOG(0,"STATUS: Calculating dipoles for current time step ...");
   std::map<atomId_t, Vector> ZeroVector;
-  BOOST_FOREACH(molecule *_mol, molecules) {
-    const Vector Dipole =
-        getDipole(
-            const_cast<const molecule *>(_mol)->begin(),
-            const_cast<const molecule *>(_mol)->end());
-    for(molecule::const_iterator iter = const_cast<const molecule *>(_mol)->begin();
-        iter != const_cast<const molecule *>(_mol)->end();
-        ++iter)
+  BOOST_FOREACH(const molecule *_mol, molecules) {
+    const Vector Dipole = getDipole(_mol->begin(),_mol->end());
+    for(molecule::const_iterator iter = _mol->begin(); iter != _mol->end(); ++iter)
       ZeroVector[(*iter)->getId()] = Dipole;
     LOG(2,"INFO: Zero alignment for molecule " << _mol->getId() << " is " << Dipole);
@@ -205,5 +200,6 @@
   World::getInstance().clearMoleculeSelection();
   World::getInstance().selectAllMolecules(MoleculeByFormula(DipoleFormula));
-  std::vector<molecule *> molecules = World::getInstance().getSelectedMolecules();
+  std::vector<const molecule *> molecules =
+      const_cast<const World &>(World::getInstance()).getSelectedMolecules();
   LOG(1,"INFO: There are " << molecules.size() << " molecules for time step " << timestep << ".");
 
@@ -212,23 +208,20 @@
   size_t i=0;
   size_t Counter_rejections = 0;
-  BOOST_FOREACH(molecule *_mol, molecules) {
-    const Vector Dipole =
-        getDipole(
-            const_cast<const molecule *>(_mol)->begin(),
-            const_cast<const molecule *>(_mol)->end());
+  BOOST_FOREACH(const molecule *_mol, molecules) {
+    const Vector Dipole = getDipole(_mol->begin(),_mol->end());
     LOG(3,"INFO: Dipole vector at time step " << timestep << " for for molecule "
         << _mol->getId() << " is " << Dipole);
     // check that all atoms are valid (zeroVector known)
-    molecule::const_iterator iter = const_cast<const molecule *>(_mol)->begin();
-    for(; iter != const_cast<const molecule *>(_mol)->end(); ++iter) {
+    molecule::const_iterator iter = _mol->begin();
+    for(; iter != _mol->end(); ++iter) {
       if (!ZeroVector.count((*iter)->getId()))
         break;
     }
-    if (iter != const_cast<const molecule *>(_mol)->end()) {
+    if (iter != _mol->end()) {
       ELOG(2, "Skipping molecule " << _mol->getName() << " as not all atoms have a valid zeroVector.");
       ++Counter_rejections;
       continue;
     } else
-      iter = const_cast<const molecule *>(_mol)->begin();
+      iter = _mol->begin();
     std::map<atomId_t, Vector>::const_iterator zeroValue = ZeroVector.find((*iter)->getId()); //due to iter is const
     double angle = 0.;
@@ -274,5 +267,6 @@
  * \return Map of doubles with values the pair of the two atoms.
  */
-DipoleCorrelationMap *DipoleCorrelation(std::vector<molecule *> &molecules)
+DipoleCorrelationMap *DipoleCorrelation(
+    const std::vector<const molecule *> &molecules)
 {
   Info FunctionInfo(__func__);
@@ -286,19 +280,13 @@
   }
 
-  for (std::vector<molecule *>::const_iterator MolWalker = molecules.begin();
+  for (std::vector<const molecule *>::const_iterator MolWalker = molecules.begin();
       MolWalker != molecules.end(); ++MolWalker) {
     LOG(2, "INFO: Current molecule is " << (*MolWalker)->getId() << ".");
-    const Vector Dipole =
-        getDipole(
-            const_cast<const molecule *>(*MolWalker)->begin(),
-            const_cast<const molecule *>(*MolWalker)->end());
-    std::vector<molecule *>::const_iterator MolOtherWalker = MolWalker;
-    for (++MolOtherWalker;
-        MolOtherWalker != molecules.end();
+    const Vector Dipole = getDipole((*MolWalker)->begin(), (*MolWalker)->end());
+    std::vector<const molecule *>::const_iterator MolOtherWalker = MolWalker;
+    for (++MolOtherWalker; MolOtherWalker != molecules.end();
         ++MolOtherWalker) {
       LOG(2, "INFO: Current other molecule is " << (*MolOtherWalker)->getId() << ".");
-      const Vector OtherDipole = getDipole(
-          const_cast<const molecule *>(*MolOtherWalker)->begin(),
-          const_cast<const molecule *>(*MolOtherWalker)->end());
+      const Vector OtherDipole = getDipole((*MolOtherWalker)->begin(), (*MolOtherWalker)->end());
       const double angle = Dipole.Angle(OtherDipole) * (180./M_PI);
       LOG(1, "Angle is " << angle << ".");
@@ -409,5 +397,8 @@
  * \return Map of dobules with values as pairs of atom and the vector
  */
-CorrelationToPointMap *CorrelationToPoint(std::vector<molecule *> &molecules, const std::vector<const element *> &elements, const Vector *point )
+CorrelationToPointMap *CorrelationToPoint(
+    const std::vector<const molecule *> &molecules,
+    const std::vector<const element *> &elements,
+    const Vector *point )
 {
   Info FunctionInfo(__func__);
@@ -421,9 +412,9 @@
   }
 
-  for (std::vector<molecule *>::const_iterator MolWalker = molecules.begin(); MolWalker != molecules.end(); MolWalker++) {
+  for (std::vector<const molecule *>::const_iterator MolWalker = molecules.begin();
+      MolWalker != molecules.end();
+      MolWalker++) {
     LOG(2, "Current molecule is " << *MolWalker << ".");
-    for (molecule::const_iterator iter = const_cast<const molecule *>(*MolWalker)->begin();
-        iter != const_cast<const molecule *>(*MolWalker)->end();
-        ++iter) {
+    for (molecule::const_iterator iter = (*MolWalker)->begin(); iter != (*MolWalker)->end(); ++iter) {
       LOG(3, "Current atom is " << **iter << ".");
       for (vector<const element *>::const_iterator type = elements.begin(); type != elements.end(); ++type)
@@ -453,5 +444,9 @@
  * \return Map of dobules with values as pairs of atom and the vector
  */
-CorrelationToPointMap *PeriodicCorrelationToPoint(std::vector<molecule *> &molecules, const std::vector<const element *> &elements, const Vector *point, const int ranges[NDIM] )
+CorrelationToPointMap *PeriodicCorrelationToPoint(
+    const std::vector<const molecule *> &molecules,
+    const std::vector<const element *> &elements,
+    const Vector *point,
+    const int ranges[NDIM] )
 {
   Info FunctionInfo(__func__);
@@ -467,11 +462,11 @@
   }
 
-  for (std::vector<molecule *>::const_iterator MolWalker = molecules.begin(); MolWalker != molecules.end(); MolWalker++) {
+  for (std::vector<const molecule *>::const_iterator MolWalker = molecules.begin();
+      MolWalker != molecules.end();
+      MolWalker++) {
     RealSpaceMatrix FullMatrix = World::getInstance().getDomain().getM();
     RealSpaceMatrix FullInverseMatrix = World::getInstance().getDomain().getMinv();
     LOG(2, "Current molecule is " << *MolWalker << ".");
-    for (molecule::const_iterator iter = const_cast<const molecule *>(*MolWalker)->begin();
-        iter != const_cast<const molecule *>(*MolWalker)->end();
-        ++iter) {
+    for (molecule::const_iterator iter = (*MolWalker)->begin(); iter != (*MolWalker)->end(); ++iter) {
       LOG(3, "Current atom is " << **iter << ".");
       for (vector<const element *>::const_iterator type = elements.begin(); type != elements.end(); ++type)
@@ -509,5 +504,9 @@
  * \return Map of doubles with values as pairs of atom and the BoundaryTriangleSet that's closest
  */
-CorrelationToSurfaceMap *CorrelationToSurface(std::vector<molecule *> &molecules, const std::vector<const element *> &elements, const Tesselation * const Surface, const LinkedCell_deprecated *LC )
+CorrelationToSurfaceMap *CorrelationToSurface(
+    const std::vector<const molecule *> &molecules,
+    const std::vector<const element *> &elements,
+    const Tesselation * const Surface,
+    const LinkedCell_deprecated *LC )
 {
   Info FunctionInfo(__func__);
@@ -522,11 +521,11 @@
   }
 
-  for (std::vector<molecule *>::const_iterator MolWalker = molecules.begin(); MolWalker != molecules.end(); MolWalker++) {
+  for (std::vector<const molecule *>::const_iterator MolWalker = molecules.begin();
+      MolWalker != molecules.end();
+      MolWalker++) {
     LOG(2, "Current molecule is " << (*MolWalker)->name << ".");
     if ((*MolWalker)->empty())
       LOG(2, "\t is empty.");
-    for (molecule::const_iterator iter = const_cast<const molecule *>(*MolWalker)->begin();
-        iter != const_cast<const molecule *>(*MolWalker)->end();
-        ++iter) {
+    for (molecule::const_iterator iter = (*MolWalker)->begin(); iter != (*MolWalker)->end(); ++iter) {
       LOG(3, "\tCurrent atom is " << *(*iter) << ".");
       for (vector<const element *>::const_iterator type = elements.begin(); type != elements.end(); ++type)
@@ -563,5 +562,10 @@
  * \return Map of doubles with values as pairs of atom and the BoundaryTriangleSet that's closest
  */
-CorrelationToSurfaceMap *PeriodicCorrelationToSurface(std::vector<molecule *> &molecules, const std::vector<const element *> &elements, const Tesselation * const Surface, const LinkedCell_deprecated *LC, const int ranges[NDIM] )
+CorrelationToSurfaceMap *PeriodicCorrelationToSurface(
+    const std::vector<const molecule *> &molecules,
+    const std::vector<const element *> &elements,
+    const Tesselation * const Surface,
+    const LinkedCell_deprecated *LC,
+    const int ranges[NDIM] )
 {
   Info FunctionInfo(__func__);
@@ -581,11 +585,11 @@
   double ShortestDistance = 0.;
   BoundaryTriangleSet *ShortestTriangle = NULL;
-  for (std::vector<molecule *>::const_iterator MolWalker = molecules.begin(); MolWalker != molecules.end(); MolWalker++) {
+  for (std::vector<const molecule *>::const_iterator MolWalker = molecules.begin();
+      MolWalker != molecules.end();
+      MolWalker++) {
     RealSpaceMatrix FullMatrix = World::getInstance().getDomain().getM();
     RealSpaceMatrix FullInverseMatrix = World::getInstance().getDomain().getMinv();
     LOG(2, "Current molecule is " << *MolWalker << ".");
-    for (molecule::const_iterator iter = const_cast<const molecule *>(*MolWalker)->begin();
-        iter != const_cast<const molecule *>(*MolWalker)->end();
-        ++iter) {
+    for (molecule::const_iterator iter = (*MolWalker)->begin(); iter != (*MolWalker)->end(); ++iter) {
       LOG(3, "Current atom is " << **iter << ".");
       for (vector<const element *>::const_iterator type = elements.begin(); type != elements.end(); ++type)
Index: src/Analysis/analysis_correlation.hpp
===================================================================
--- src/Analysis/analysis_correlation.hpp	(revision a58c16f47065e9446b7916a407271c6544180022)
+++ src/Analysis/analysis_correlation.hpp	(revision 99db9baca1782e208fb17d2cc6ace3f4cadb568b)
@@ -61,30 +61,73 @@
 /********************************************** declarations *******************************/
 
-range<size_t> getMaximumTrajectoryBounds(const std::vector<atom *> &atoms);
-std::map<atomId_t, Vector> CalculateZeroAngularDipole(const std::vector<molecule *> &molecules);
-
-DipoleAngularCorrelationMap *DipoleAngularCorrelation(const Formula &DipoleFormula, const size_t timestep, const std::map<atomId_t, Vector> &ZeroVector, const enum ResetWorldTime DoTimeReset = DontResetTime);
-DipoleCorrelationMap *DipoleCorrelation(std::vector<molecule *> &molecules);
+range<size_t> getMaximumTrajectoryBounds(
+    const std::vector<const atom *> &atoms);
+std::map<atomId_t, Vector> CalculateZeroAngularDipole(
+    const std::vector<const molecule *> &molecules);
+
+DipoleAngularCorrelationMap *DipoleAngularCorrelation(
+    const Formula &DipoleFormula,
+    const size_t timestep,
+    const std::map<atomId_t, Vector> &ZeroVector,
+    const enum ResetWorldTime DoTimeReset = DontResetTime);
+DipoleCorrelationMap *DipoleCorrelation(
+    const std::vector<const molecule *> &molecules);
 PairCorrelationMap *PairCorrelation(
     const World::ConstAtomComposite &atoms_first,
     const World::ConstAtomComposite &atoms_second,
     const double max_distance);
-CorrelationToPointMap *CorrelationToPoint(std::vector<molecule *> &molecules, const std::vector<const element *> &elements, const Vector *point );
-CorrelationToSurfaceMap *CorrelationToSurface(std::vector<molecule *> &molecules, const std::vector<const element *> &elements, const Tesselation * const Surface, const LinkedCell_deprecated *LC );
-CorrelationToPointMap *PeriodicCorrelationToPoint(std::vector<molecule *> &molecules, const std::vector<const element *> &elements, const Vector *point, const int ranges[NDIM] );
-CorrelationToSurfaceMap *PeriodicCorrelationToSurface(std::vector<molecule *> &molecules, const std::vector<const element *> &elements, const Tesselation * const Surface, const LinkedCell_deprecated *LC, const int ranges[NDIM] );
-int GetBin ( const double value, const double BinWidth, const double BinStart );
-void OutputCorrelation_Header( ofstream * const file );
-void OutputCorrelation_Value( ofstream * const file, BinPairMap::const_iterator &runner );
-void OutputDipoleAngularCorrelation_Header( ofstream * const file );
-void OutputDipoleAngularCorrelation_Value( ofstream * const file, DipoleAngularCorrelationMap::const_iterator &runner );
-void OutputDipoleCorrelation_Header( ofstream * const file );
-void OutputDipoleCorrelation_Value( ofstream * const file, DipoleCorrelationMap::const_iterator &runner );
-void OutputPairCorrelation_Header( ofstream * const file );
-void OutputPairCorrelation_Value( ofstream * const file, PairCorrelationMap::const_iterator &runner );
-void OutputCorrelationToPoint_Header( ofstream * const file );
-void OutputCorrelationToPoint_Value( ofstream * const file, CorrelationToPointMap::const_iterator &runner );
-void OutputCorrelationToSurface_Header( ofstream * const file );
-void OutputCorrelationToSurface_Value( ofstream * const file, CorrelationToSurfaceMap::const_iterator &runner );
+CorrelationToPointMap *CorrelationToPoint(
+    const std::vector<const molecule *> &molecules,
+    const std::vector<const element *> &elements,
+    const Vector *point );
+CorrelationToSurfaceMap *CorrelationToSurface(
+    const std::vector<const molecule *> &molecules,
+    const std::vector<const element *> &elements,
+    const Tesselation * const Surface,
+    const LinkedCell_deprecated *LC );
+CorrelationToPointMap *PeriodicCorrelationToPoint(
+    const std::vector<const molecule *> &molecules,
+    const std::vector<const element *> &elements,
+    const Vector *point, const int ranges[NDIM] );
+CorrelationToSurfaceMap *PeriodicCorrelationToSurface(
+    const std::vector<const molecule *> &molecules,
+    const std::vector<const element *> &elements,
+    const Tesselation * const Surface,
+    const LinkedCell_deprecated *LC,
+    const int ranges[NDIM] );
+int GetBin (
+    const double value,
+    const double BinWidth,
+    const double BinStart );
+void OutputCorrelation_Header(
+    ofstream * const file );
+void OutputCorrelation_Value(
+    ofstream * const file,
+    BinPairMap::const_iterator &runner );
+void OutputDipoleAngularCorrelation_Header(
+    ofstream * const file );
+void OutputDipoleAngularCorrelation_Value(
+    ofstream * const file,
+    DipoleAngularCorrelationMap::const_iterator &runner );
+void OutputDipoleCorrelation_Header(
+    ofstream * const file );
+void OutputDipoleCorrelation_Value(
+    ofstream * const file,
+    DipoleCorrelationMap::const_iterator &runner );
+void OutputPairCorrelation_Header(
+    ofstream * const file );
+void OutputPairCorrelation_Value(
+    ofstream * const file,
+    PairCorrelationMap::const_iterator &runner );
+void OutputCorrelationToPoint_Header(
+    ofstream * const file );
+void OutputCorrelationToPoint_Value(
+    ofstream * const file,
+    CorrelationToPointMap::const_iterator &runner );
+void OutputCorrelationToSurface_Header(
+    ofstream * const file );
+void OutputCorrelationToSurface_Value(
+    ofstream * const file,
+    CorrelationToSurfaceMap::const_iterator &runner );
 
 /** Searches for lowest and highest value in a given map of doubles.
Index: src/Analysis/unittests/AnalysisCorrelationToPointUnitTest.cpp
===================================================================
--- src/Analysis/unittests/AnalysisCorrelationToPointUnitTest.cpp	(revision a58c16f47065e9446b7916a407271c6544180022)
+++ src/Analysis/unittests/AnalysisCorrelationToPointUnitTest.cpp	(revision 99db9baca1782e208fb17d2cc6ace3f4cadb568b)
@@ -105,5 +105,5 @@
   // init maps
   World::getInstance().selectAllMolecules(AllMolecules());
-  allMolecules = World::getInstance().getSelectedMolecules();
+  allMolecules = const_cast<const World &>(World::getInstance()).getSelectedMolecules();
   CPPUNIT_ASSERT_EQUAL( (size_t) 1, allMolecules.size());
   pointmap = CorrelationToPoint( allMolecules, elements, (const Vector *)point );
Index: src/Analysis/unittests/AnalysisCorrelationToPointUnitTest.hpp
===================================================================
--- src/Analysis/unittests/AnalysisCorrelationToPointUnitTest.hpp	(revision a58c16f47065e9446b7916a407271c6544180022)
+++ src/Analysis/unittests/AnalysisCorrelationToPointUnitTest.hpp	(revision 99db9baca1782e208fb17d2cc6ace3f4cadb568b)
@@ -41,5 +41,5 @@
 private:
 
-      std::vector<molecule *> allMolecules;
+      std::vector<const molecule *> allMolecules;
       molecule *TestMolecule;
       const element *hydrogen;
Index: src/Analysis/unittests/AnalysisCorrelationToSurfaceUnitTest.cpp
===================================================================
--- src/Analysis/unittests/AnalysisCorrelationToSurfaceUnitTest.cpp	(revision a58c16f47065e9446b7916a407271c6544180022)
+++ src/Analysis/unittests/AnalysisCorrelationToSurfaceUnitTest.cpp	(revision 99db9baca1782e208fb17d2cc6ace3f4cadb568b)
@@ -145,5 +145,5 @@
 
   World::getInstance().selectAllMolecules(AllMolecules());
-  allMolecules = World::getInstance().getSelectedMolecules();
+  allMolecules = const_cast<const World &>(World::getInstance()).getSelectedMolecules();
   CPPUNIT_ASSERT_EQUAL( (size_t) 2, allMolecules.size());
 
Index: src/Analysis/unittests/AnalysisCorrelationToSurfaceUnitTest.hpp
===================================================================
--- src/Analysis/unittests/AnalysisCorrelationToSurfaceUnitTest.hpp	(revision a58c16f47065e9446b7916a407271c6544180022)
+++ src/Analysis/unittests/AnalysisCorrelationToSurfaceUnitTest.hpp	(revision 99db9baca1782e208fb17d2cc6ace3f4cadb568b)
@@ -49,5 +49,5 @@
 private:
 
-      std::vector<molecule *> allMolecules;
+      std::vector<const molecule *> allMolecules;
       molecule *TestSurfaceMolecule;
       const element *hydrogen;
Index: src/Parser/TremoloParser.cpp
===================================================================
--- src/Parser/TremoloParser.cpp	(revision a58c16f47065e9446b7916a407271c6544180022)
+++ src/Parser/TremoloParser.cpp	(revision 99db9baca1782e208fb17d2cc6ace3f4cadb568b)
@@ -727,9 +727,12 @@
 }
 
-bool FormatParser< tremolo >::saveAtomsInExttypes(std::ostream &output, const std::vector<atom*> &atoms, const int id) const
+bool FormatParser< tremolo >::saveAtomsInExttypes(
+    std::ostream &output,
+    const std::vector<const atom*> &atoms,
+    const int id) const
 {
   bool status = true;
   // parse the file
-  for (std::vector<atom *>::const_iterator iter = atoms.begin();
+  for (std::vector<const atom *>::const_iterator iter = atoms.begin();
       iter != atoms.end(); ++iter) {
     const int atomicid = getLocalId((*iter)->getId());
Index: src/Parser/TremoloParser.hpp
===================================================================
--- src/Parser/TremoloParser.hpp	(revision a58c16f47065e9446b7916a407271c6544180022)
+++ src/Parser/TremoloParser.hpp	(revision 99db9baca1782e208fb17d2cc6ace3f4cadb568b)
@@ -85,5 +85,8 @@
    * @return true - all atomic ids are valid, false - one id was not set (-1)
    */
-  bool saveAtomsInExttypes(std::ostream &output, const std::vector<atom*> &atoms, const int id) const;
+  bool saveAtomsInExttypes(
+      std::ostream &output,
+      const std::vector<const atom*> &atoms,
+      const int id) const;
 
 protected:
Index: src/Shapes/ShapeFactory.cpp
===================================================================
--- src/Shapes/ShapeFactory.cpp	(revision a58c16f47065e9446b7916a407271c6544180022)
+++ src/Shapes/ShapeFactory.cpp	(revision 99db9baca1782e208fb17d2cc6ace3f4cadb568b)
@@ -138,5 +138,6 @@
     case MoleculeSurfaceType:
     {
-      const std::vector<molecule *> molecules = World::getInstance().getSelectedMolecules();
+      const std::vector<const molecule *> molecules =
+          const_cast<const World &>(World::getInstance()).getSelectedMolecules();
       if (molecules.empty()) {
         ELOG(2, "No molecule is selected, cannot factorize its bounding shape.");
Index: src/World.cpp
===================================================================
--- src/World.cpp	(revision a58c16f47065e9446b7916a407271c6544180022)
+++ src/World.cpp	(revision 99db9baca1782e208fb17d2cc6ace3f4cadb568b)
@@ -701,10 +701,21 @@
 }
 
-const std::vector<atom *> World::getSelectedAtoms() const {
+std::vector<atom *> World::getSelectedAtoms() {
   std::vector<atom *> returnAtoms;
-  returnAtoms.resize(countSelectedAtoms());
-  int count = 0;
-  for (AtomSet::const_iterator iter = selectedAtoms.begin(); iter != selectedAtoms.end(); ++iter)
-    returnAtoms[count++] = iter->second;
+  std::transform(
+      selectedAtoms.begin(),
+      selectedAtoms.end(),
+      back_inserter(returnAtoms),
+      _take<atom*,World::AtomSet::value_type>::get);
+  return returnAtoms;
+}
+
+std::vector<const atom *> World::getSelectedAtoms() const {
+  std::vector<const atom *> returnAtoms;
+  std::transform(
+      selectedAtoms.begin(),
+      selectedAtoms.end(),
+      back_inserter(returnAtoms),
+      _take<atom*,World::AtomSet::value_type>::get);
   return returnAtoms;
 }
Index: src/World.hpp
===================================================================
--- src/World.hpp	(revision a58c16f47065e9446b7916a407271c6544180022)
+++ src/World.hpp	(revision 99db9baca1782e208fb17d2cc6ace3f4cadb568b)
@@ -429,5 +429,6 @@
   bool isSelected(const atom *_atom) const;
   bool isAtomSelected(const atomId_t no) const;
-  const std::vector<atom *> getSelectedAtoms() const;
+  std::vector<atom *> getSelectedAtoms();
+  std::vector<const atom *> getSelectedAtoms() const;
 
   void clearMoleculeSelection();
