/*
    Copyright (C) 2012 Bartosz Broda, Jacek Bochenek, Łukasz Biłeńkij, 
    Damian Jaworski, Wojciech Mazur, Dominik Piasecki, Maciej Piasecki, 
    Radosław Ramocki
    
    This program is free software; you can redistribute it and/or modify it
under the terms of the GNU Lesser 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 LICENSE and COPYING files for more details.
*/

#include "ROCK.h"

namespace clustering
{
	namespace ROCK
	{
		ROCK::ROCK(boost::shared_ptr<Matrices::SuperMatrix> pMatrix,
			boost::shared_ptr<smartcomparator::simfunctions::LocalSimilarityFunction> simfun,
			double simlinktreshold)
			: m_pMatrix(pMatrix),
			m_pSimilarityFunction(simfun),
			m_dSimLinkTreshold(simlinktreshold),
			m_bFlatClustering(false),
			m_pOutputClusters(new ClusteringBase::ClusterSet()),
			m_stNumberOfGroups(5),
			m_stMinSizeOfCluster(1),
			m_pClusterSet(new ClusteringBase::ClusterSet()),
			m_pClusterMap(new ClusteringBase::ClusterMap())
//			m_pOutputLinks(new ClusteringBase::LinkSet()),
		{

		}

		ClusteringBase::ClusterSetPtr ROCK::GetClusters()
		{
			return m_pOutputClusters;
		}

/*		ClusteringBase::LinkSetPtr ROCK::GetLinks()
		{
			return m_pOutputLinks;
		}*/

		void ROCK::CountAlgorithmParams()
		{
			m_dFi = (1 - m_dSimLinkTreshold) / (1 + m_dSimLinkTreshold);
			m_dFi = (1+m_dFi)/2;

			std::cout << "m_dFi = " << m_dFi << std::endl;
		}

		void ROCK::DrawRandomSample()
		{
			std::cout << "Drawing random sample..." << std::endl;
			for(size_t i = 0; i < m_pMatrix->NumberOfRows(); i++)
			{
				m_vSampleIndices.push_back(i);
			}
			std::cout << "Drawing random sample...Done!" << std::endl;
		}//void DrawRandomSample()


		void ROCK::StartClustering()
		{
			DrawRandomSample();
			CountAlgorithmParams();
			// utworzenie grup 1-el

			CreateSingeltonClusters();

			CreateInitialLinks();
#ifdef CLUSTERING_DEBUG_ROCK
			TestLinksDifferences();
#endif



			std::cout << "Clustering..." << std::endl;
			Cluster();
			LabelReminingData();

		}//void ROCK::StartClustering()

#ifdef CLUSTERING_DEBUG_ROCK
		void ROCK::TestLinksDifferences()
		{
			using namespace ClusteringBase;
			size_t diffNo = 0;
			for(ClusterSet::iterator i1 = m_pClusterSet->begin(); i1 != m_pClusterSet->end(); i1++)
			{
				size_t ind1 = (*i1)->GetID();
				for(ClusterSet::iterator i2 = m_pClusterSet->begin(); i2 != m_pClusterSet->end(); i2++)
				{
					size_t ind2 = (*i2)->GetID();
					if (GetNumberOfLinks(ind1, ind2) != GetNumberOfLinks(ind2, ind1))
					{
						diffNo++;
					}
				}
			}

			if(diffNo)
				std::cerr << "Liczba roznych: " << diffNo << std::endl;
		}
#endif

		void ROCK::Cluster()
		{
			using namespace ClusteringBase;

			std::cout << "Creating local heaps...\n";
			// utworzenie lokalnych kopcow
			CreateLocalHeaps();
			std::cout << "Creating local heaps... Done!\n";

			// utworznie globalnego kopca
			BuildGlobalHeap();

			// wlasciwe grupowanie 
			size_t iterNo = 0;

			while(m_GlobalHeap.Size()> m_stNumberOfGroups)
			{
				//std::cout << "Global heap size: " << m_GlobalHeap.Size() << std::endl;
#ifdef CLUSTERING_DEBUG_ROCK
				TestLinksDifferences();
#endif
				iterNo++;
				// wyciagniecie max 
				ClusterWithGoodness cwgU = m_GlobalHeap.PopMaxElement();//6
				ClusterPtr u = cwgU.pCluster;//6
				size_t u_id = u->GetID();

				// wyciagniecie najlepszego dla u
				ClusterWithGoodness cwgV = m_LocalHeaps[u_id].GetMaxElement();

				// usuniecie v (czyli nie ma juz na stosie globalnym ani u ani v)
				m_GlobalHeap.RemoveElement(cwgV.pCluster);
				ClusterPtr v = cwgV.pCluster;//7
				size_t v_id = v->GetID();

				//scalenie i dodanei grupy
				ClusterPtr w = MergeClusters( u,v) ;//9
				AddCluster(w);//dodany z ifa
/*				if(!m_bFlatClustering)
				{
					
					AddCluster(w);
					if(u->GetPointsID().size() >= m_stMinSizeOfCluster &&
						v->GetPointsID().size() >= m_stMinSizeOfCluster)
					{
						LinkPtr l1(new Link(u_id,w->GetID() ));
						LinkPtr l2(new Link(v_id,w->GetID()));
						m_pOutputLinks->insert(l1);
						m_pOutputLinks->insert(l2);
					}

				}*/

				std::set<size_t> indices = m_LocalHeaps[u_id].GetAllClusterIndices();
				std::set<size_t> &vIndices = m_LocalHeaps[v_id].GetAllClusterIndices();
				indices.insert(vIndices.begin(), vIndices.end());

				for (std::set<size_t>::iterator i = indices.begin();
					i!= indices.end(); i++)
				{
					if(*i == u_id || *i == v_id)
						continue;

					ClusterMap::iterator find = m_pClusterMap->find(*i);
					if(find == m_pClusterMap->end())
						continue;
					AddMergedCluster(w, find->second, u,v);
				}

				ClusterWithGoodness cwgW = m_LocalHeaps[w->GetID()].GetMaxElement();

				if(m_LocalHeaps[w->GetID()].Size() > 0)
				{
					//std::cout << "------------------- DODAJE\n";
					m_GlobalHeap.AddOrUpdateElement(cwgW.pCluster, cwgW.dGoodness);
				}

				m_LocalHeaps.erase(u_id);//17
				m_LocalHeaps.erase(v_id);//17
				//m_Links.erase(u_id);
				//m_Links.erase(v_id);
				if(iterNo % 10 == 0)
					std::cout << "ROCK Iteration: " << iterNo << std::endl;
			}

			std::cout << "Wielkosc kopca Q = " << m_GlobalHeap.Size() << std::endl;
/*			if(m_bFlatClustering)
			{
				m_pOutputLinks->clear();
				m_pOutputClusters->clear();

				while (m_GlobalHeap.Size()> 0)
				{
					ClusterWithGoodness cwgU = m_GlobalHeap.PopMaxElement();
					ClusterPtr u = cwgU.pCluster;
					AddCluster(u);
				}
			}*/

		}//void Cluster();

		void ROCK::AddMergedCluster(ClusteringBase::ClusterPtr w, ClusteringBase::ClusterPtr x, 
			ClusteringBase::ClusterPtr u, ClusteringBase::ClusterPtr v)
		{

			using namespace ClusteringBase;
			size_t u_id = u->GetID();
			size_t v_id = v->GetID();
			size_t x_id = x->GetID();
			size_t w_id = w->GetID();


			int iLinkXU = GetNumberOfLinks(x_id, u_id);
			int iLinkXV = GetNumberOfLinks(x_id, v_id);


			if (iLinkXU != GetNumberOfLinks(u_id, x_id))
				std::cerr << "**** **** **** niesymetryczne linki...\n";

			if (iLinkXV != GetNumberOfLinks(v_id, x_id))
				std::cerr << "* * * * * * niesymetryczne linki...\n";

			int sum = iLinkXU + iLinkXV;
			if(sum > 0)
			{
				//std::map<size_t, int> &xLinks = m_Links[x_id];
				//std::map<size_t, int> &wLinks = m_Links[w_id];

				m_Links[x_id][w_id] = sum;
				m_Links[w_id][x_id] = sum;


				//xLinks[w_id] = sum;//11
				//wLinks[x_id] = sum;//11
			}

			ClusterWithGoodness oldBestForX = m_LocalHeaps[x_id].GetMaxElement();

			m_LocalHeaps[x_id].RemoveElement(u);//12
			m_LocalHeaps[x_id].RemoveElement(v);//12

			double g_x_w = CountGoodness(x,w);

			if(g_x_w > 0)
			{
				m_LocalHeaps[x_id].AddOrUpdateElement(w,g_x_w);//13
				m_LocalHeaps[w_id].AddOrUpdateElement(x,g_x_w);//13
			}

			//XXX: bylo &newBestForX -> sprawdzic jaki to moze miec wplyw
			ClusterWithGoodness newBestForX  = m_LocalHeaps[x_id].GetMaxElement();

			//14:

			if(newBestForX.pCluster.get() == NULL && oldBestForX.pCluster.get() == NULL)
				return;

			if(newBestForX.pCluster.get() == NULL && oldBestForX.pCluster.get() != NULL)
			{
				m_GlobalHeap.RemoveElement(oldBestForX.pCluster);
				return;
			}

			if(oldBestForX.pCluster.get() == NULL && newBestForX.pCluster.get() != NULL)
			{
				m_GlobalHeap.AddOrUpdateElement(newBestForX.pCluster, newBestForX.dGoodness);
				return;
			}

			if(oldBestForX.pCluster->GetID() != newBestForX.pCluster->GetID())
			{
				//m_GlobalHeap.RemoveElement(oldBestForX.pCluster);
				m_GlobalHeap.AddOrUpdateElement(newBestForX.pCluster, newBestForX.dGoodness);
			}//*/
		}


		void ROCK::CreateLocalHeaps()
		{
			for(ClusteringBase::ClusterSet::iterator c = m_pClusterSet->begin(); c != m_pClusterSet->end();
				c++)
				BuildLocalHeap(*c);

		}//	void ROCK::CreateLocalHeaps()

		void ROCK::BuildLocalHeap(ClusteringBase::ClusterPtr cluster_i)
		{
			using namespace ClusteringBase;

			size_t id = cluster_i->GetID();

			for(ClusterSet::iterator cluster_j = m_pClusterSet->begin(); 
				cluster_j != m_pClusterSet->end(); cluster_j++)
			{
				if(id == (*cluster_j)->GetID())
					continue;

				double goodness = CountGoodness(cluster_i, *cluster_j);
				if(goodness < 0.000001)
					continue;

				//ClusterWithGoodness cwg(*cluster_j, goodness);

				m_LocalHeaps[id].AddOrUpdateElement(*cluster_j, goodness);

			}


		}//BuildLocalHeap

		void ROCK::BuildGlobalHeap()
		{

			using namespace ClusteringBase;
			for (ClusterSet::iterator cluster = m_pClusterSet->begin(); 
				cluster != m_pClusterSet->end(); cluster++)
			{
				ClusterPtr pCluster = (*cluster);

				if (m_LocalHeaps[pCluster ->GetID()].Size() < 1)
					continue;

				ClusterWithGoodness cwg = m_LocalHeaps[pCluster ->GetID()].GetMaxElement();

				double bestGoodnes = cwg.dGoodness;

				m_GlobalHeap.AddOrUpdateElement(pCluster, bestGoodnes);

			}//*/
		}//void ROCK::BuildGlobalHeap()


		void ROCK::CreateSingeltonClusters()
		{
			std::cout << "Creating singelton clusters..." << std::endl;
			for(size_t sampleindex = 0; sampleindex < m_vSampleIndices.size(); sampleindex++)
			{
				size_t docindex = m_vSampleIndices[sampleindex];
				ClusteringBase::ClusterPtr cluster(new ClusteringBase::Cluster() );
				cluster->AddPoint(docindex);
				AddCluster(cluster);
				m_InitialClusterToSample[cluster->GetID()] = sampleindex;
			}//for each sample index
			std::cout << "Creating singelton clusters...Done!" << std::endl;

		}//void ROCK::CreateSingeltonClusters()

		void ROCK::CountNeighbors(std::vector< std::vector<size_t> > &n)
		{
			using namespace ClusteringBase;
			using Matrices::CompressedVector;
			std::cout << "Counting initial neighbors...\n";
			double avg_size = 0;

			size_t size = m_pClusterSet->size();
			size_t i = 1;

			for (ClusterSet::iterator clusterFromIter = m_pClusterSet->begin();
				clusterFromIter != m_pClusterSet->end(); clusterFromIter++, i++)
			{
				ClusterPtr clusterFrom = *clusterFromIter;
				size_t fromSampleIndex = m_InitialClusterToSample[clusterFrom->GetID()];
				Matrices::CompressedVectorPtr  fromDoc = m_pMatrix->GetCRS()->GetCompressedRow(fromSampleIndex);

				std::vector<size_t> fromN;
				n.push_back(fromN);
				std::vector<size_t> &fromNeighbors = n[n.size()-1];

				for (ClusterSet::iterator clusterToIter = m_pClusterSet->begin(); 
					clusterToIter != m_pClusterSet->end(); clusterToIter++)
				{
					ClusterPtr clusterTo = *clusterToIter;
					size_t toSampleIndex = m_InitialClusterToSample[clusterTo->GetID()];
					if(fromSampleIndex == toSampleIndex)
						continue;

					Matrices::CompressedVectorPtr toDoc = m_pMatrix->GetCRS()->GetCompressedRow(toSampleIndex);
					double sim = m_pSimilarityFunction->CountSimilarity(*fromDoc, *toDoc);

					sim = (sim+1)/2;

					if(sim >= m_dSimLinkTreshold)
						fromNeighbors.push_back(clusterTo->GetID());
					
				}
				avg_size+= fromNeighbors.size();
				if( i % 10 == 0)
					std::cout << "Counting neighbors, " << i << " / " << size << std::endl;

			}


			avg_size /= size;

			std::cout << "\nCounting initial neighbors... DONE\n";
			std::cout << "Avarage number of neighbours: " << avg_size << std::endl;
		}//CountNeighbors

		void ROCK::CreateInitialLinks()
		{
			std::vector< std::vector<size_t> > neighbors;
			CountNeighbors(neighbors);
			m_InitialClusterToSample.clear();

			std::cout << "Counting initial links..." << std::endl;
			size_t size = neighbors.size();
			for(size_t i= 0; i < size ; i++)
			{
				std::vector<size_t> &N = neighbors[i];

				if(N.size() == 0)
					continue;

				for(size_t j = 0; j < N.size() -1; j++)
				{
					for(size_t l =j+1; l < N.size(); l++)
					{
						// jako ze powiazania sa pomiedzy grupami, a sasiedztwo pomiedzy
						// dokumentami, to trzeba zamienic...
						size_t N_j = N [ j ];
						size_t N_l = N [ l ];

						m_Links[ N_j ][N_l]++;
						m_Links[ N_l ][N_j]++;
					}
				}
				if( i % 100 == 0 )
					std::cout << "Creating initial links, " << i << " / " << size << std::endl;

			}
			std::cout << "\nCounting initial links... Done" << std::endl;

		}//void CreateLinkMatrix();


		int ROCK::GetNumberOfLinks(size_t c1, size_t c2)
		{
			std::map<size_t, std::map<size_t, int> >::iterator find = m_Links.find(c1);
			if(find == m_Links.end())
				return 0;

			std::map<size_t, int>::iterator findValue = find->second.find(c2);

			if(findValue == find->second.end())
				return 0;

			return findValue->second;

		}
		double ROCK::CountGoodness(ClusteringBase::ClusterPtr  C1, 
			ClusteringBase::ClusterPtr C2)
		{
			using namespace ClusteringBase;


			double link_Ci_Cj = (double)GetNumberOfLinks(C1->GetID(), C2->GetID());

			if(link_Ci_Cj == 0)
				return 0;

			double n_i  = (double) C1->GetPointsID().size();
			double n_j  = (double) C2->GetPointsID().size();

			double pow_factor = 1. + 2.*m_dFi;

			double goodness = link_Ci_Cj / 
				( pow(n_i+n_j, pow_factor) - pow(n_i, pow_factor) - pow(n_j, pow_factor) );

			return goodness;
		}//double ROCK::CountGoodness


		void ROCK::AddCluster(ClusteringBase::ClusterPtr &pCluster)
		{
			m_pClusterSet->insert(pCluster);
			m_pClusterMap->insert( ClusteringBase::IdClusterPtrPair
				(pCluster->GetID(), pCluster)	);
			if(pCluster->GetPointsID().size() >= m_stMinSizeOfCluster)
				m_pOutputClusters->insert(pCluster);
		}//void ROCK::AddCluster(ClusteringBase::ClusterPtr &pCluster)


		void ROCK::LabelReminingData()
		{
		}//void LabelDataOnDisk();

		void ROCK::CreateParentsMap()
		{
		    std::cout << "Creating parents map" << std::endl;
		    for (size_t i = 0; i < m_vSampleIndices.size(); ++i)
			std::cout << m_vSampleIndices[i] << std::endl;

		    ClusteringBase::ClusterSet leafs;
		    ClusteringBase::GetLeafs(this->GetClusters(), leafs);

		    if (leafs.size() != m_pMatrix->NumberOfRows())
			throw std::runtime_error("Number of leafs not equals to number of rows"); //should be thrown if the random sample is used

		    for (ClusteringBase::ClusterSet::iterator iter = leafs.begin(); iter != leafs.end(); ++iter)
			this->FindParent(*iter);

		    std::cout << "Done..." << std::endl;
		}

		void ROCK::FindParent(ClusteringBase::ClusterPtr cluster)
		{
		    if (cluster->GetParent() != NULL)
		    {
			m_ParentsMap[cluster->GetID()] = cluster->GetParent()->GetID();
			FindParent(cluster->GetParent());
		    }
		    else
		    {
			m_ParentsMap[cluster->GetID()] = -1;
			return;
		    }
		}

		void ROCK::SaveClutoTreeFile(const std::string & path)
		{
		    std::ofstream file( path.c_str() ) ;
		    if( !file )
			    throw std::runtime_error( "Error while opening file: " + path );

		    for (std::map<size_t, int>::const_iterator iter = m_ParentsMap.begin(); iter != m_ParentsMap.end(); ++iter)
		    {
			file << iter->second << " -1 -1" << std::endl;
			std::cout << iter->first << " " << iter->second << std::endl;
		    }
		    file.close();
		}

	}//namespace ROCK



}//namespace Clustering
