/*
 * Project: MoleCuilder
 * Description: creates and alters molecular systems
 * Copyright (C)  2012 University of Bonn. All rights reserved.
 * Please see the COPYING file or "Copyright notice" in builder.cpp for details.
 *
 *
 *   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/>.
 */

/*
 * Extractors.cpp
 *
 *  Created on: 15.10.2012
 *      Author: heber
 */

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

#include "CodePatterns/MemDebug.hpp"

#include <utility>
#include <vector>
#include <boost/assign.hpp>

#include "CodePatterns/Log.hpp"

#include "LinearAlgebra/Vector.hpp"

#include "FunctionApproximation/Extractors.hpp"
#include "FunctionApproximation/FunctionArgument.hpp"

using namespace boost::assign;

FunctionModel::arguments_t
Extractors::_detail::gatherAllDistanceArguments(
    const Fragment::positions_t &positions,
    const size_t globalid)
{
  FunctionModel::arguments_t result;

  // go through current configuration and gather all other distances
  Fragment::positions_t::const_iterator firstpositer = positions.begin();
  for (;firstpositer != positions.end(); ++firstpositer) {
    Fragment::positions_t::const_iterator secondpositer = positions.begin();//firstpositer;
    for (; secondpositer != positions.end(); ++secondpositer) {
      if (firstpositer == secondpositer)
        continue;
      argument_t arg;
      const Vector firsttemp((*firstpositer)[0],(*firstpositer)[1],(*firstpositer)[2]);
      const Vector secondtemp((*secondpositer)[0],(*secondpositer)[1],(*secondpositer)[2]);
      arg.distance = firsttemp.distance(secondtemp);
      arg.indices = std::make_pair(
          std::distance(
              positions.begin(), firstpositer),
          std::distance(
              positions.begin(), secondpositer)
          );
      arg.globalid = globalid;
      result.push_back(arg);
    }
  }

  return result;
}

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::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);
}