/*    This file is part of the Gudhi Library. The Gudhi library
 *    (Geometric Understanding in Higher Dimensions) is a generic C++
 *    library for computational topology.
 *
 *    Author(s):       David Salinas
 *
 *    Copyright (C) 2014  INRIA Sophia Antipolis-Mediterranee (France)
 *
 *    This program 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 3 of the License, or
 *    (at your option) any later version.
 *
 *    This program 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 this program.  If not, see <http://www.gnu.org/licenses/>.
 */

#ifndef GUDHI_SKELETON_BLOCKER_CONTRACTOR_H_
#define GUDHI_SKELETON_BLOCKER_CONTRACTOR_H_

#include <memory>
#include <cassert>

// todo remove the queue to be independent from cgald
#include "gudhi/Contraction/CGAL_queue/Modifiable_priority_queue.h"
//#include <CGAL/Modifiable_priority_queue.h>


#include <list>
#include <boost/scoped_array.hpp>
#include <boost/scoped_ptr.hpp>


#include "gudhi/Contraction/Edge_profile.h"
#include "gudhi/Contraction/policies/Cost_policy.h"
#include "gudhi/Contraction/policies/Edge_length_cost.h"
#include "gudhi/Contraction/policies/Placement_policy.h"
#include "gudhi/Contraction/policies/First_vertex_placement.h"

#include "gudhi/Contraction/policies/Valid_contraction_policy.h"
#include "gudhi/Contraction/policies/Dummy_valid_contraction.h" //xxx remove
#include "gudhi/Contraction/policies/Link_condition_valid_contraction.h" //xxx remove

#include "gudhi/Contraction/policies/Contraction_visitor.h"

#include "gudhi/Skeleton_blocker/Skeleton_blocker_complex_visitor.h"

#include "gudhi/Utils.h"

namespace Gudhi{

namespace contraction {


template <class Profile>
Placement_policy<Profile>* make_first_vertex_placement(){
	return new First_vertex_placement<Profile>();
}

template <class Profile>
Valid_contraction_policy<Profile>* make_link_valid_contraction(){
	return new Link_condition_valid_contraction<Profile>();
}


// visitor that remove popable blockers after an edge contraction
template <class Profile>
class Contraction_visitor_remove_popable : public Contraction_visitor<Profile>{
public:
	typedef typename Profile::Point Point;
	typedef typename Profile::Complex Complex;
	typedef typename Complex::Vertex_handle Vertex_handle;

	//todo explore only neighborhood with stack
	void on_contracted(const Profile  &profile, boost::optional< Point > placement) override{
		std::list<Vertex_handle> vertex_to_check;
		vertex_to_check.push_front(profile.v0_handle());

		while(!vertex_to_check.empty()){
			Vertex_handle v = vertex_to_check.front();
			vertex_to_check.pop_front();

			bool blocker_popable_found=true;
			while (blocker_popable_found){
				blocker_popable_found = false;
				for(auto block : profile.complex().blocker_range(v)){
					if (profile.complex().is_popable_blocker(block)) {
						for(Vertex_handle nv : *block)
							if(nv!=v) vertex_to_check.push_back(nv);
						profile.complex().delete_blocker(block);
						blocker_popable_found = true;
						break;
					}
				}
			}
		}
	}
};

template <class Profile>
Contraction_visitor<Profile>* make_remove_popable_blockers_visitor(){
	return new Contraction_visitor_remove_popable<Profile>();
}






/**
 *@class Skeleton_blocker_contractor
 *@brief Class that allows to contract iteratively edges of a simplicial complex.
 *
 * @details Basically, the simplification algorithm consists in iteratively picking the
 * edge with lowest cost and performing an edge contraction if the contraction is valid.
 * This class is policy based (and much inspired from the edge collapse package of CGAL).
 *
 * Policies that can be changed are :
 *  - the cost policy : how much cost an edge contraction
 *  - the placement policy : where will be placed the contraction point
 *  - the valid contraction policy : is the contraction valid. For instance, it can be
 *  a topological condition (link condition) or a geometrical condition (normals messed up).
 *
 * TODO expliquer la pile
 */
template<
class GeometricSimplifiableComplex,
class EdgeProfile = Edge_profile<GeometricSimplifiableComplex>
>
class Skeleton_blocker_contractor : private skbl::Dummy_complex_visitor<typename GeometricSimplifiableComplex::Vertex_handle>{

	GeometricSimplifiableComplex& complex_;

public:
	typedef typename GeometricSimplifiableComplex::Graph_vertex Graph_vertex;
	typedef typename GeometricSimplifiableComplex::Vertex_handle Vertex_handle;
	typedef typename GeometricSimplifiableComplex::Simplex_handle Simplex_handle;
	typedef typename GeometricSimplifiableComplex::Simplex_handle_iterator Simplex_handle_iterator;



	typedef typename GeometricSimplifiableComplex::Root_vertex_handle Root_vertex_handle;

	typedef typename GeometricSimplifiableComplex::Graph_edge Graph_edge;
	typedef typename GeometricSimplifiableComplex::Edge_handle Edge_handle;
	typedef typename GeometricSimplifiableComplex::Point Point;

	typedef EdgeProfile  Profile;


	typedef Cost_policy<Profile> Cost_policy_;
	typedef Placement_policy<Profile> Placement_policy_;
	typedef Valid_contraction_policy<Profile> Valid_contraction_policy_;
	typedef Contraction_visitor<EdgeProfile> Contraction_visitor_;
	typedef Edge_profile_factory<EdgeProfile> Edge_profile_factory_;



	typedef boost::optional<double> Cost_type;
	typedef boost::optional<Point> Placement_type ;

	typedef size_t size_type;

	typedef Skeleton_blocker_contractor Self ;



private:

	struct Compare_id
	{
		Compare_id() : algorithm_(0) {}

		Compare_id( Self const* aAlgorithm ) : algorithm_(aAlgorithm) {}

		bool operator() ( Edge_handle a, Edge_handle b ) const
		{
			return algorithm_->get_undirected_edge_id(a) < algorithm_->get_undirected_edge_id(b);
		}

		Self const* algorithm_ ;
	} ;

	struct Compare_cost
	{
		Compare_cost() : algorithm_(0) {
		}

		Compare_cost( Self const* aAlgorithm ) : algorithm_(aAlgorithm) {
		}

		bool operator() ( Edge_handle a, Edge_handle b ) const
		{
			// NOTE: A cost is an optional<> value.
			// Absent optionals are ordered first; that is, "none < T" and "T > none" for any defined T != none.
			// In consequence, edges with undefined costs will be promoted to the top of the priority queue and poped out first.
			return algorithm_->get_data(a).cost() < algorithm_->get_data(b).cost();
		}

		Self const* algorithm_ ;

	} ;

	struct Undirected_edge_id : boost::put_get_helper<size_type, Undirected_edge_id>
	{
		typedef boost::readable_property_map_tag category;
		typedef size_type                        value_type;
		typedef size_type                        reference;
		typedef Edge_handle                  key_type;

		Undirected_edge_id() : algorithm_(0) {}

		Undirected_edge_id( Self const* aAlgorithm ) : algorithm_(aAlgorithm) {}

		size_type operator[] ( Edge_handle e ) const { return algorithm_->get_undirected_edge_id(e); }

		Self const* algorithm_ ;
	} ;

	typedef CGAL::Modifiable_priority_queue<Edge_handle,Compare_cost,Undirected_edge_id> PQ ;
	typedef typename PQ::handle pq_handle ;


	// An Edge_data is associated with EVERY edge in the complex (collapsible or not).
	// It relates the edge with the PQ-handle needed to update the priority queue
	// It also relates the edge with a policy-based cache
	class Edge_data
	{
	public :

		Edge_data() : PQHandle_(),cost_() {}

		Cost_type const& cost() const { return cost_ ; }
		Cost_type      & cost()       { return cost_ ; }

		pq_handle PQ_handle() const { return PQHandle_ ;}

		bool is_in_PQ() const { return PQHandle_ != PQ::null_handle() ; }

		void set_PQ_handle( pq_handle h ) { PQHandle_ = h ; }

		void reset_PQ_handle() { PQHandle_ = PQ::null_handle() ; }

	private:
		pq_handle PQHandle_ ;
		Cost_type cost_ ;

	} ;
	typedef Edge_data* Edge_data_ptr ;
	typedef boost::scoped_array<Edge_data> Edge_data_array ;


	int get_undirected_edge_id ( Edge_handle edge ) const {
		return complex_[edge].index() ;
	}


	const Edge_data& get_data ( Edge_handle edge ) const
	{
		return edge_data_array_[get_undirected_edge_id(edge)];
	}


	Edge_data& get_data ( Edge_handle edge )
	{
		return edge_data_array_[get_undirected_edge_id(edge)];
	}

	Cost_type get_cost(const Profile & profile) const{
		return (*cost_policy_)(profile,get_placement(profile));
	}

	Profile create_profile(Edge_handle edge) const{
		if(edge_profile_factory_)
			return edge_profile_factory_->make_profile(complex_,edge);
		else
			return Profile(complex_,edge);
	}


	void insert_in_PQ( Edge_handle edge, Edge_data& data )
	{
		data.set_PQ_handle(heap_PQ_->push(edge));
		++current_num_edges_heap_;
	}

	void update_in_PQ( Edge_handle edge, Edge_data& data )
	{
		data.set_PQ_handle(heap_PQ_->update(edge,data.PQ_handle())) ;
	}

	void remove_from_PQ( Edge_handle edge, Edge_data& data )
	{
		data.set_PQ_handle(heap_PQ_->erase(edge,data.PQ_handle()));
		--current_num_edges_heap_;
	}

	boost::optional<Edge_handle> pop_from_PQ()	{
		boost::optional<Edge_handle> edge = heap_PQ_->extract_top();
		if ( edge )
		{
			get_data(*edge).reset_PQ_handle();
		}
		return edge ;
	}



private:

	/**
	 * @brief Collect edges.
	 *
	 * Iterates over all edges of the simplicial complex and
	 * 1) inserts them in the priority queue sorted according to the Cost policy.
	 * 2) set the id() field of each edge
	 */
	void collect_edges(){
		//
		// Loop over all the edges in the complex in the heap
		//
		size_type size = complex_.num_edges();
		DBG("Collecting edges ...");
		DBGMSG("num edges :",size);

		edge_data_array_.reset( new Edge_data[size] ) ;

		heap_PQ_.reset( new PQ (size, Compare_cost(this), Undirected_edge_id(this) ) ) ;

		std::size_t id = 0 ;

		//xxx do a parralel for
		for(auto edge : complex_.edge_range()){
			//			Edge_handle edge = *edge_it;
			complex_[edge].index() = id++;
			Profile const& profile = create_profile(edge);
			Edge_data& data = get_data(edge);
			data.cost() = get_cost(profile) ;
			++initial_num_edges_heap_;
			insert_in_PQ(edge,data);
			if(contraction_visitor_) contraction_visitor_->on_collected(profile,data.cost());
		}

		DBG("Edges collected.");

	}

	bool should_stop(double lCost,const Profile &profile) const{
		return false;
	}

	boost::optional<Point> get_placement(const Profile& profile) const{
		return (*placement_policy_)(profile);
	}

	bool is_contraction_valid( Profile const& profile, Placement_type placement ) const{
		return (*valid_contraction_policy_)(profile,placement);
	}


public:
	/**
	 * \brief Contract edges.
	 *
	 * While the heap is not empty, it extracts the edge with the minimum
	 * cost in the heap then try to contract it.
	 * It stops when the Stop policy says so or when the number of contractions
	 * given by 'num_max_contractions' is reached (if this number is positive).
	 */
	void contract_edges(int num_max_contractions=-1){

		DBG("\n\nContract edges");
		int num_contraction = 0 ;

		bool unspecified_num_contractions = (num_max_contractions == -1);
		//
		// Pops and processes each edge from the PQ
		//
		boost::optional<Edge_handle> edge ;
		while ( (edge = pop_from_PQ())&& ((num_contraction<num_max_contractions)||(unspecified_num_contractions)))
		{
			Profile const& profile = create_profile(*edge);
			Cost_type cost(get_data(*edge).cost());
			if(contraction_visitor_) contraction_visitor_->on_selected(profile,cost,0,0);

			DBGMSG("\n\n---- Pop edge - num vertices :",complex_.num_vertices());

			if (cost){
				DBGMSG("sqrt(cost):",std::sqrt(*cost));
				if (should_stop(*cost,profile) )
				{
					if(contraction_visitor_) contraction_visitor_->on_stop_condition_reached(profile);
					DBG("should_stop");
					break ;
				}
				Placement_type placement = get_placement(profile);
				if ( is_contraction_valid(profile,placement) && placement )
				{
					DBG("contraction_valid");
					contract_edge(profile,placement);
					++ num_contraction;
				}
				else
				{
					DBG("contraction not valid");
					if(contraction_visitor_) contraction_visitor_->on_non_valid(profile);
				}
			}
			else
			{
				DBG("uncomputable cost");
			}
		}
	}


	bool is_in_heap(Edge_handle edge) const{
		if(heap_PQ_->empty()) return false;
		else{
			return edge_data_array_[get_undirected_edge_id(edge)].is_in_PQ();
		}
	}

	bool is_heap_empty() const{
		return heap_PQ_->empty();
	}


	/**
	 * @brief Returns an Edge_handle and a Placement_type. This pair consists in
	 * the edge with the lowest cost in the heap together with its placement.
	 * The returned value is initialized iff the heap is non-empty.
	 */
	boost::optional<std::pair<Edge_handle,Placement_type > > top_edge(){
		boost::optional<std::pair<Edge_handle,Placement_type > > res;

		if(!heap_PQ_->empty()) {
			auto edge = heap_PQ_->top();
			Profile const& profile = create_profile(edge);
			Placement_type placement = get_placement(profile);
			res = make_pair(edge,placement);
			DBGMSG("top edge:",complex_[edge]);

		}
		return res;
	}


	/**
	 * @brief Constructor with default policies.
	 *
	 * @details The default cost, placement, valid and visitor policies
	 * are respectively : the edge length, the first point, the link condition
	 */
	Skeleton_blocker_contractor(GeometricSimplifiableComplex& complex)
	:complex_(complex),
	 cost_policy_(new Edge_length_cost<Profile>),
	 placement_policy_(new First_vertex_placement<Profile>),
	 valid_contraction_policy_(new Link_condition_valid_contraction<Profile>),
	 contraction_visitor_(new Contraction_visitor_()),
	 edge_profile_factory_(0),
	 initial_num_edges_heap_(0),
	 current_num_edges_heap_(0)
	{
		complex_.set_visitor(this);
		if(contraction_visitor_) contraction_visitor_->on_started(complex);
		collect_edges();
	}

	/**
	 * @brief Constructor with customed policies.
	 * @remark Policies destruction is handle by the class with smart pointers.
	 */
	Skeleton_blocker_contractor(GeometricSimplifiableComplex& complex,
			Cost_policy_ *cost_policy,
			Placement_policy_ * placement_policy = new First_vertex_placement<Profile>,
			Valid_contraction_policy_ * valid_contraction_policy = new Link_condition_valid_contraction<Profile>,
			Contraction_visitor_* contraction_visitor = new Contraction_visitor_(),
			Edge_profile_factory_* edge_profile_factory = NULL
	):
		complex_(complex),
		cost_policy_(cost_policy),
		placement_policy_(placement_policy),
		valid_contraction_policy_(valid_contraction_policy),
		contraction_visitor_(contraction_visitor),
		edge_profile_factory_(edge_profile_factory),
		initial_num_edges_heap_(0),
		current_num_edges_heap_(0)
	{
		complex_.set_visitor(this);
		if(contraction_visitor) contraction_visitor->on_started(complex);
		collect_edges();
	}


	~Skeleton_blocker_contractor(){
		complex_.set_visitor(0);
	}

private:


	void contract_edge(const Profile& profile, Placement_type placement ) {
		if(contraction_visitor_) contraction_visitor_->on_contracting(profile,placement);

		assert(complex_.contains_vertex(profile.v0_handle()));
		assert(complex_.contains_vertex(profile.v1_handle()));
		assert(placement);

		profile.complex().point(profile.v0_handle()) = *placement;

		// remark : this is not necessary since v1 will be deactivated
		//	profile.complex().point(profile.v1_handle()) = *placement;

		complex_.contract_edge(profile.v0_handle(),profile.v1_handle());

		assert(complex_.contains_vertex(profile.v0_handle()));
		assert(!complex_.contains_vertex(profile.v1_handle()));

		update_changed_edges();

		// the visitor could do something as complex_.remove_popable_blockers();
		if(contraction_visitor_) contraction_visitor_->on_contracted(profile,placement);
	}

private:

	// every time the visitor's method on_changed_edge is called, it adds an
	// edge to changed_edges_
	std::vector< Edge_handle > changed_edges_;

	/**
	 * @brief we update the cost and the position in the heap of an edge that has
	 * been changed
	 */
	inline void on_changed_edge(Vertex_handle a,Vertex_handle b) override{
		boost::optional<Edge_handle> ab(complex_[std::make_pair(a,b)]);
		assert(ab);
		changed_edges_.push_back(*ab);
	}

	void update_changed_edges(){
		//xxx do a parralel for

		DBG("update edges");

		//		sequential loop
		for(auto ab : changed_edges_){
			//1-get the Edge_handle corresponding to ab
			//2-change the data in mEdgeArray[ab.id()]
			//3-update the heap
			Edge_data& data = get_data(ab);
			Profile const& profile = create_profile(ab);
			data.cost() = get_cost(profile) ;
			if ( data.is_in_PQ()){
				update_in_PQ(ab,data);
			}
			else{
				insert_in_PQ(ab,data);
			}
		}
		changed_edges_.clear();
	}


private:
	void on_remove_edge(Vertex_handle a,Vertex_handle b) override{

		boost::optional<Edge_handle> ab((complex_[std::make_pair(a,b)]));
		assert(ab);
		Edge_data& lData = get_data(*ab) ;
		if ( lData.is_in_PQ() )
		{
			remove_from_PQ(*ab,lData) ;
		}
	}
private:
	/**
	 * @brief Called when the edge 'ax' has been added while the edge 'bx'
	 * is still there but will be removed on next instruction.
	 * We assign the index of 'bx' to the edge index of 'ax'
	 */
	void on_swaped_edge(Vertex_handle a,Vertex_handle b,Vertex_handle x) override{
		boost::optional<Edge_handle> ax(complex_[std::make_pair(a,x)]);
		boost::optional<Edge_handle> bx(complex_[std::make_pair(b,x)]);
		assert(ax&& bx);
		complex_[*ax].index() =complex_[*bx].index();
	}
private:
	/**
	 * @brief Called when a blocker is removed.
	 * All the edges that passes through the blocker may be edge-contractible
	 * again and are thus reinserted in the heap.
	 */
	void on_delete_blocker(const Simplex_handle * blocker) override{
		// we go for all pairs xy that belongs to the blocker
		// note that such pairs xy are necessarily edges of the complex
		// by definition of a blocker

		// todo uniqument utile pour la link condition
		// laisser à l'utilisateur? booleen update_heap_on_removed_blocker?
		Simplex_handle blocker_copy(*blocker);
		for (auto x = blocker_copy.begin(); x!= blocker_copy.end(); ++x){
			for(auto y=x ; ++y != blocker_copy.end(); ){
				auto edge_descr(complex_[std::make_pair(*x,*y)]);
				assert(edge_descr);
				Edge_data& data = get_data(*edge_descr);
				Profile const& profile = create_profile(*edge_descr);
				data.cost() = get_cost(profile) ;

				// If the edge is already in the heap
				// its priority has not changed.
				// If the edge is not present, we reinsert it
				// remark : we could also reinsert the edge
				// only if it is valid
				if ( !data.is_in_PQ() ){
					insert_in_PQ(*edge_descr,data);
				}
			}
		}
	}


private:
	std::shared_ptr<Cost_policy_> cost_policy_;
	std::shared_ptr<Placement_policy_> placement_policy_;
	std::shared_ptr<Valid_contraction_policy_> valid_contraction_policy_;
	std::shared_ptr<Contraction_visitor_> contraction_visitor_;

	//in case the user wants to do something special when the edge profile
	//are created (for instance add some info)
	std::shared_ptr<Edge_profile_factory_> edge_profile_factory_;
	Edge_data_array edge_data_array_ ;

	boost::scoped_ptr<PQ> heap_PQ_ ;
	int initial_num_edges_heap_;
	int current_num_edges_heap_;

};

}  // namespace contraction
} // namespace GUDHI

#endif /* GUDHI_SKELETON_BLOCKER_CONTRACTOR_H_ */