Index: src/FunctionApproximation/Extractors.cpp
===================================================================
--- src/FunctionApproximation/Extractors.cpp	(revision 4a77be7c74807798597478c1e2a88d09178390fe)
+++ src/FunctionApproximation/Extractors.cpp	(revision 301dbfe7ce5eed79d11e438d8f8b33e350ce1480)
@@ -40,4 +40,5 @@
 #include <boost/assign.hpp>
 
+#include "CodePatterns/Assert.hpp"
 #include "CodePatterns/Log.hpp"
 
@@ -81,43 +82,120 @@
 }
 
-FunctionModel::arguments_t
-Extractors::gatherFirstDistance(
-    const Fragment& fragment,
-    const size_t index,
-    const size_t firstelement,
-    const size_t secondelement
-    ) {
-  const Fragment::charges_t charges = fragment.getCharges();
-  const Fragment::positions_t positions = fragment.getPositions();
-  typedef Fragment::charges_t::const_iterator chargeiter_t;
-  std::vector< chargeiter_t > firstpair;
-  firstpair.reserve(2);
-  firstpair +=
-      std::find(charges.begin(), charges.end(), firstelement),
-      std::find(charges.begin(), charges.end(), secondelement);
-  if ((firstpair[0] == charges.end()) || (firstpair[1] == charges.end())) {
-    // complain if tuple not found
-    ELOG(1, "Could not find pair " << firstelement << "," << secondelement
-        << " in fragment " << fragment);
-    return FunctionModel::arguments_t();
-  }
-  // convert position_t to Vector
-  std::vector< std::pair<Vector, size_t> > DistancePair;
-  for (std::vector<chargeiter_t>::const_iterator firstpairiter = firstpair.begin();
-      firstpairiter != firstpair.end(); ++firstpairiter) {
+Fragment::positions_t Extractors::_detail::gatherPositionsFromCharges(
+    const Fragment::positions_t &positions,
+    const Fragment::charges_t &charges,
+    const chargeiters_t targets
+    )
+{
+  Fragment::positions_t filtered_positions;
+  for (chargeiters_t::const_iterator firstpairiter = targets.begin();
+      firstpairiter != targets.end(); ++firstpairiter) {
     Fragment::positions_t::const_iterator positer = positions.begin();
     const size_t steps = std::distance(charges.begin(), *firstpairiter);
     std::advance(positer, steps);
-    DistancePair.push_back(
-        std::make_pair(Vector((*positer)[0], (*positer)[1], (*positer)[2]),
-            steps));
-  }
-  // finally convert Vector pair to distance-like argument
-  argument_t arg;
-  arg.indices.first = DistancePair[0].second;
-  arg.indices.second = DistancePair[1].second;
-  arg.distance = DistancePair[0].first.distance(DistancePair[1].first);
-  arg.globalid = index;
-
-  return FunctionModel::arguments_t(1, arg);
-}
+    filtered_positions.push_back(*positer);
+  }
+  return filtered_positions;
+}
+
+Fragment::positions_t
+Extractors::gatherDistanceOfTuples(
+    const Fragment& fragment,
+    const Fragment::charges_t elements
+    ) {
+  const Fragment::charges_t charges = fragment.getCharges();
+  /// The main problem here is that we have to know how many same
+  /// elements (but different atoms!) we are required to find. Hence,
+  /// we first have to count same elements, then get different targets
+  /// for each and then associated them in correct order back again.
+
+  // 1. we have to make elements unique with counts, hence convert to map
+  elementcounts_t elementcounts;
+  for (Fragment::charges_t::const_iterator elementiter = elements.begin();
+      elementiter != elements.end(); ++elementiter) {
+    // insert new element
+    std::pair< elementcounts_t::iterator, bool> inserter =
+        elementcounts.insert( std::make_pair( *elementiter, 1) );
+    // if already present, just increase its count
+    if (!inserter.second)
+      ++(inserter.first->second);
+  }
+
+  // 2. then for each element we need as many targets (chargeiters) as counts
+  elementtargets_t elementtargets;
+  for (elementcounts_t::const_iterator countiter = elementcounts.begin();
+      countiter != elementcounts.end();
+      ++countiter) {
+    chargeiter_t chargeiter = charges.begin();
+    const element_t &element = countiter->first;
+    const count_t &count = countiter->second;
+    for (count_t i = 0; i < count; ++i) {
+      chargeiter_t tempiter = std::find(chargeiter, charges.end(), element);
+      if (tempiter != charges.end()) {
+        // try to insert new list
+        std::pair< elementtargets_t::iterator, bool> inserter =
+          elementtargets.insert( std::make_pair( countiter->first, chargeiters_t(1, tempiter)) );
+        // if already present, append to it
+        if (!inserter.second) {
+          inserter.first->second.push_back(tempiter);
+        } else { // if created, increase vector's reserve to known size
+          inserter.first->second.reserve(countiter->second);
+        }
+        // search from this element onwards then
+        chargeiter = ++tempiter;
+      } else {
+        ELOG(1, "Could not find desired number of elements " << elements << " in fragment.");
+        return Fragment::positions_t();
+      }
+    }
+  }
+
+  // 3. we go again through elements and use one found target for each count
+  // in that order
+  elementcounts_t counts; // how many chargeiters of this element have been used
+  chargeiters_t targets;
+  targets.reserve(elements.size());
+  for (Fragment::charges_t::const_iterator elementiter = elements.begin();
+      elementiter != elements.end(); ++elementiter) {
+    const element_t &element = *elementiter;
+    count_t &count = counts[element]; // if not present, std::map creates instances with default of 0
+#ifndef NDEBUG
+    {
+      elementcounts_t::const_iterator testiter = elementcounts.find(element);
+      ASSERT( (testiter != elementcounts.end()) && (count < testiter->second),
+          "Extractors::gatherDistanceOfTuples() - we want to use more chargeiters for element "
+          +toString(element)+" than we counted initially.");
+    }
+#endif
+    elementtargets_t::iterator targetiter = elementtargets.find(element);
+    ASSERT (targetiter != elementtargets.end(),
+        "Extractors::gatherDistanceOfTuples() - not enough chargeiters for element "
+        +toString(element)+".");
+    chargeiters_t &chargeiters = targetiter->second;
+    const chargeiter_t &chargeiter = chargeiters[count++];
+    targets.push_back(chargeiter);
+  }
+#ifndef NDEBUG
+  // check all for debugging
+  for (chargeiters_t::const_iterator chargeiter = targets.begin();
+      chargeiter != targets.end();
+      ++chargeiter)
+    ASSERT( *chargeiter != charges.end(),
+        "Extractors::gatherDistanceOfTuples() - we have not found enough targets?!");
+#endif
+  // 4. convert position_t to Vector
+  return Extractors::_detail::gatherPositionsFromCharges(
+          fragment.getPositions(),
+          charges,
+          targets);
+}
+
+FunctionModel::arguments_t Extractors::reorderArgumentsByIncreasingDistance(
+    const FunctionModel::arguments_t &args
+    )
+{
+  FunctionModel::arguments_t returnargs(args);
+  std::sort(returnargs.begin(), returnargs.end(), argument_t::DistanceComparator);
+  return returnargs;
+}
+
Index: src/FunctionApproximation/Extractors.hpp
===================================================================
--- src/FunctionApproximation/Extractors.hpp	(revision 4a77be7c74807798597478c1e2a88d09178390fe)
+++ src/FunctionApproximation/Extractors.hpp	(revision 301dbfe7ce5eed79d11e438d8f8b33e350ce1480)
@@ -23,4 +23,12 @@
  */
 namespace Extractors {
+  typedef Fragment::charges_t::const_iterator chargeiter_t;
+  typedef std::vector<chargeiter_t> chargeiters_t;
+
+  typedef size_t count_t;
+  typedef Fragment::charge_t element_t;
+  typedef std::map< element_t, count_t> elementcounts_t;
+  typedef std::map< element_t, chargeiters_t > elementtargets_t;
+
   /** Namespace for some internal helper functions.
    *
@@ -37,4 +45,20 @@
         const Fragment::positions_t &positions,
         const size_t globalid);
+
+    /** Gather all positions from the same aligned vector of charges.
+     *
+     * Basically, we filter the positions indicated by the targets but
+     * from a different vector that has the same layout.
+     *
+     * \param positions positions to search
+     * \param charges charges to match with \targets
+     * \param targets iterators on charges
+     * \return filtered positions
+     */
+    Fragment::positions_t gatherPositionsFromCharges(
+        const Fragment::positions_t &positions,
+        const Fragment::charges_t &charges,
+        const chargeiters_t targets
+        );
   }
 
@@ -53,17 +77,22 @@
   }
 
-  /** Gather first distance for the two matching charges.
+  /** Gather all positions associated to the matching \a elements.
    *
    * \param fragment fragment with all nuclei positions
-   * \param index index refers to the index within the global set of configurations
-   * \param firstelement first element of pair
-   * \param secondelement second element of pair, order is reflected in indices of return argument_t
-   * \return vector of of argument_t containing all found distances
+   * \param elements tuple of desired elements
+   * \return vector of positions_t containing
    */
-  FunctionModel::arguments_t gatherFirstDistance(
+  Fragment::positions_t gatherDistanceOfTuples(
       const Fragment& fragment,
-      const size_t index,
-      const size_t firstelement,
-      const size_t secondelement
+      const Fragment::charges_t elements
+      );
+
+  /** Reorder arguments by increasing distance.
+   *
+   * \param args arguments to reorder
+   * \return reordered args
+   */
+  FunctionModel::arguments_t reorderArgumentsByIncreasingDistance(
+      const FunctionModel::arguments_t &args
       );
 
Index: src/LevMartester.cpp
===================================================================
--- src/LevMartester.cpp	(revision 4a77be7c74807798597478c1e2a88d09178390fe)
+++ src/LevMartester.cpp	(revision 301dbfe7ce5eed79d11e438d8f8b33e350ce1480)
@@ -364,5 +364,9 @@
     // Afterwards we go through all of this type and gather the distance and the energy value
     TrainingData MorseData(
-        boost::bind(&Extractors::gatherFirstDistance, _1, _2, 6, 6) // gather first carbon pair
+        boost::bind(&Extractors::_detail::gatherAllDistanceArguments,
+            boost::bind(&Extractors::gatherDistanceOfTuples,
+                _1, Fragment::charges_t(2,6.)
+            ), _2 // gather first carbon pair
+          )
         );
     MorseData(homologies.getHomologousGraphs(graph));
