/*
    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 "cbc.h"
#include <iostream>
#include "kbest.h"
#include "clutowrapper.h"
#include "config.h"

#include "boost/date_time/posix_time/posix_time.hpp"

void kbest_test();

int main(int argc, char* argv[])
{
	using namespace clustering::cbc;

	Options opt(argc, argv);

	switch (Config::remove_features)
	{
	case NORMAL: std::cout << "Normal removing features...\n"; break;
	case MIN: std::cout << "Remove min mode...\n"; break;
	case WEIGHTED: std::cout << "weighted removing features..\n"; break;
	}
	if (Config::verbose)
		std::cout << "Verbose mode ON" << std::endl;
	CBC cbc(false); // true = only phase 3
	cbc.StartCBC();
//	std::string f = "/home/boban/committees_i2_rwf";
//	cbc.ResumePhase3( f );
	//kbest_test();
	return 0;
}

void kbest_test()
{
	using namespace clustering::cbc;
	std::vector<std::string>  words;
	words.push_back("zamek");
	words.push_back("agent");
	words.push_back("zgoda");
	words.push_back("dom");
	{
		// following are for generating SimMtrix
		Matrices::SuperMatrixPtr pMatrix = Matrices::IO::LoadSuperMatrix(Config::dir,Config::matrix_name, false, true);
        smartcomparator::weights::RankBuilderPtr pRWF(new smartcomparator::weights::NormalRankBuilder(pMatrix) );
		smartcomparator::compare::CompareFunctionPtr pCmp ( new smartcomparator::compare::RWFZScoreRankComparision(pMatrix,pRWF) );
		KBest::GenerateFromMatrix( Config::kbest_file, 1000, 0.0001, pMatrix, pCmp );
		exit(0);
	}


	KBestPtr pKBest = KBestPtr(new KBest() );
	Matrices::SymmetricSimMatrixPtr pSimMatrix = pKBest->ReadFromGeneratedFiles( Config::kbest_file );
	Matrices::SuperMatrixPtr pMatrix = Matrices::IO::LoadSuperMatrix(Config::dir,Config::matrix_name, false, true);
    smartcomparator::weights::RankBuilderPtr pRWF(new smartcomparator::weights::NormalRankBuilder(pMatrix) );
	smartcomparator::compare::CompareFunctionPtr pCmp ( new smartcomparator::compare::RWFZScoreRankComparision(pMatrix, pRWF) );
	smartcomparator::compare::CompareFunctionPtr pCmpMI ( new smartcomparator::compare::PMIComparision(pMatrix, false) );

	smartcomparator::compare::Comparator comparator(pMatrix,pCmp);
	smartcomparator::compare::Comparator comparator_mi(pMatrix,pCmpMI);
	for(size_t w = 0; w < words.size(); w++)
	{
		std::string &word = words[w];
		std::vector<std::string> &kbest = pKBest->GetKBest(word);
		std::cout << word << "\n";
		std::cout << "\tMI, from file...\n";
		for(size_t i = 0; i < 20 && i < kbest.size(); i++)
			std::cout << "\t\t" << kbest[i] << std::endl;

		std::cout << "\n\tMI without discounting, from matrix...\n";
		std::multimap<double, size_t>  best_mi = comparator_mi.FindSimiliar( pMatrix->GetRowIndex( word ), 20 );
		for(std::multimap<double, size_t>::iterator iter = best_mi.begin(); iter != best_mi.end(); iter++)
			std::cout <<"\t\t"<< pMatrix->GetRowFeatures()->operator[]( iter->second ).GetLabel() << std::endl;

		std::cout << "\n\tZScore, from matrix...\n";
		std::multimap<double, size_t>  best = comparator.FindSimiliar( pMatrix->GetRowIndex( word ), 20 );
		for(std::multimap<double, size_t>::iterator iter = best.begin(); iter != best.end(); iter++)
			std::cout <<"\t\t"<< pMatrix->GetRowFeatures()->operator[]( iter->second ).GetLabel() << std::endl;
	}

}

namespace clustering
{
	namespace cbc
	{
		CBC::CBC(bool only_phase_3)
		{
			m_pMatrix = Matrices::IO::LoadSuperMatrix(Config::dir,Config::matrix_name, false, true);
			smartcomparator::weights::RankBuilderPtr pRWF(new smartcomparator::weights::NormalRankBuilder(m_pMatrix) );
			m_pCmp = smartcomparator::compare::CompareFunctionPtr( new smartcomparator::compare::RWFZScoreRankComparision(m_pMatrix,pRWF) );

	//Added reading kbest - test it!
			KBest::GenerateFromMatrix(Config::kbest_file, 1000, 0.0001, m_pMatrix, m_pCmp);
			if(!only_phase_3)
			{
				m_pKBest = KBestPtr(new KBest() );
				m_pSimMatrix = m_pKBest->ReadFromGeneratedFiles( Config::kbest_file, "", Config::simmatrix_top_x );

				size_t N = 0;
				for(size_t row_no = 0; row_no < m_pMatrix->NumberOfRows(); row_no++)
					N += (* (m_pMatrix->GetRowFeatures()))[row_no].GetTermFrequency();
				m_dN = (double) (N);
			}
			m_bUseTScoreRank = true;
			if(m_bUseTScoreRank)
			{
                smartcomparator::weights::RankBuilderPtr pRWF(new smartcomparator::weights::NormalRankBuilder(m_pMatrix) );
				boost::shared_ptr<smartcomparator::weights::RWFZScore > pWeight(new
								  smartcomparator::weights::RWFZScore(m_pMatrix, pRWF) );
				smartcomparator::weights::MatrixTransformer mt(pWeight);
				mt.Transform();
			}
			else
			{
				boost::shared_ptr<smartcomparator::weights::SimpleMIWeight > pWeight(new
								  smartcomparator::weights::SimpleMIWeight(m_pMatrix) );
				smartcomparator::weights::MatrixTransformer mt(pWeight);
				mt.Transform();
			}
			std::cout << "Commitee output file: " << Config::committee_output_file << std::endl;
		}

		void CBC::StartCBC()
		{
			Phase2();
			Phase3();
			WriteResults();
		}

		void CBC::WriteResults()
		{
			WriteWordsToCommmittees();
			WriteCommitteesToWords();
			if (Config::cluto_file != "")
				SaveClutoFile();
		}

		void CBC::WriteCommitteesToWords()
		{
			std::string filename = Config::committee_output_file + ".committees";
			std::cout << "Writing results to file " << filename << std::endl;
			std::ofstream file( filename.c_str() );
			std::ostream *out = &file;

			if(!out)
			{
				std::cerr << "Cannot open output file: " << Config::committee_output_file << "\nWriting to std::out\n";
				out = &std::cout;
			}

			for(std::map<size_t, std::set<WordSimPair, wsp_sim_cmp> >::iterator iter = m_mCommitteeIndexToWord.begin();
					iter != m_mCommitteeIndexToWord.end(); iter++)
			{
				(*out) << "( committee no: " << iter->first << " = ";
				Committee& committee  = m_Committees[iter->first];
				std::copy( committee.second.begin(), committee.second.end(), std::ostream_iterator<std::string>( *out, ", ") );
				(*out) << std::endl;

				for(std::set<WordSimPair, wsp_sim_cmp>::iterator citer = iter->second.begin(); citer != iter->second.end(); citer++)
					(*out) << "\t " << (*citer).word << "\n";
				(*out) << ")" << std::endl;
			}

			out->flush();
			std::cout << "Writing results to file " << Config::committee_output_file << "... Done! " << std::endl;
		}

		void CBC::SaveClutoFile()
		{
			std::string filename = Config::cluto_file;
			std::cout << "Writing results in CLUTO format to file " << filename << std::endl;
			std::ofstream file( filename.c_str() );
			std::ostream *out = &file;

			if(!out)
			{
				std::cerr << "Cannot open output file: " << Config::cluto_file << "\nWriting to std::out\n";
				out = &std::cout;
			}

			for (size_t i = 0; i < m_pMatrix->NumberOfRows(); ++i)
			{
				const std::string& word = m_pMatrix->GetRowFeatures()->operator[](i).GetLabel();
				std::map<std::string, std::set<size_t> >::iterator iter = m_mWordToCommitteeIndex.find(word);
				if (iter == m_mWordToCommitteeIndex.end())
				{
					(*out) << (-1) << std::endl;
					continue;
				}

				for(std::set<size_t>::iterator citer = iter->second.begin(); citer != iter->second.end(); citer++)
				{
					(*out) << *citer << " ";
					//Committee& committee  = m_Committees[*citer];
					//std::copy( committee.second.begin(), committee.second.end(), std::ostream_iterator<std::string>( *out, ", ") );
				}
				(*out) << std::endl;
			}

			out->flush();
			std::cout << "Writing results to file " << Config::cluto_file << "... Done! " << std::endl;
		}

		void CBC::WriteWordsToCommmittees()
		{
			std::string filename = Config::committee_output_file + ".words";
			std::cout << "Writing results to file " << filename << std::endl;
			std::ofstream file( filename.c_str() );
			std::ostream *out = &file;

			if(!out)
			{
				std::cerr << "Cannot open output file: " << Config::committee_output_file << "\nWriting to std::out\n";
				out = &std::cout;
			}

			for(std::map<std::string, std::set<size_t> >::iterator iter = m_mWordToCommitteeIndex.begin();
					iter != m_mWordToCommitteeIndex.end(); iter++)
			{
				if(iter->second.size() < 2 )
					continue;
				(*out) << "(" << iter->first << std::endl;
				for(std::set<size_t>::iterator citer = iter->second.begin(); citer != iter->second.end(); citer++)
				{
					(*out) << "\t " << *citer<< "\t ";
					Committee& committee  = m_Committees[*citer];
					std::copy( committee.second.begin(), committee.second.end(), std::ostream_iterator<std::string>( *out, ", ") );
					(*out) << std::endl;
				}
				(*out) << ")" << std::endl;
			}

			for(std::map<std::string, std::set<size_t> >::iterator iter = m_mWordToCommitteeIndex.begin();
					iter != m_mWordToCommitteeIndex.end(); iter++)
			{
				if(iter->second.size() > 1  )
					continue;
				(*out) << "(" << iter->first << std::endl;
				for(std::set<size_t>::iterator citer = iter->second.begin(); citer != iter->second.end(); citer++)
				{
					(*out) << "\t " << *citer<< "\t ";
					Committee& committee  = m_Committees[*citer];
					std::copy( committee.second.begin(), committee.second.end(), std::ostream_iterator<std::string>( *out, ", ") );
					(*out) << std::endl;
				}
				(*out) << ")" << std::endl;
			}

			out->flush();
			std::cout << "Writing results to file " << Config::committee_output_file << "... Done! " << std::endl;
		}

		void CBC::Phase2()
		{
			std::cout << "CBC: entring Phase2" << std::endl;
			ClutoWrapper cluto(m_pMatrix, m_pCmp, m_pKBest, m_pSimMatrix);

			//select words for clustering... (TODO: addittion of reading words from file -> extracting to different method?)
			std::vector<std::string> words;
			for(size_t rowno = 0; rowno < m_pMatrix->NumberOfRows(); rowno++)
			{
//				std::cout << "TF: " << m_pMatrix->GetRowFeatures()->operator[](rowno).GetTermFrequency() << std::endl;
//				std::cout << "NZ: " << m_pMatrix->GetRowFeatures()->operator[](rowno).GetNZ() << std::endl;
				if(m_pMatrix->GetRowFeatures()->operator[](rowno).GetTermFrequency() > Config::committee_min_tf
						&& m_pMatrix->GetRowFeatures()->operator[](rowno).GetNZ() > Config::committee_min_nz )
				{
					Matrices::CompressedVectorPtr pRow = m_pMatrix->GetRow(rowno);
					double sum_tf = Measures::CountVectorSum( pRow );
					if( sum_tf > Config::committee_sum_tf_threshold)
						words.push_back( m_pMatrix->GetRowFeatures()->operator[](rowno).GetLabel() );
				}
			}

			std::cout << "Selected " << words.size() << " words with sum_tf > " << Config::committee_sum_tf_threshold << ", minnz = " << Config::committee_min_nz
				<< " mintf = " << Config::committee_min_tf << std::endl;

			size_t iteration_number = 0;
			while(words.size() )//while we have something to do in phase 2 (recursive step)
			{
				std::cout << "ENTER ITERATION" << std::endl;
				std::copy(words.begin(), words.end(), std::ostream_iterator<std::string>(std::cout, " "));
				// step 1 - cluster top similar lists
				cluto.StartClustering(words, Config::workdir+boost::lexical_cast<std::string>(iteration_number), Config::k );
				ClustersPtr pClusters = cluto.GetClusters();
				std::cout << "Cluster size = " << pClusters->size() << std::endl;
				Committees current_committee;
				size_t current_check = 0;
				// step 2-3 select commitiee only, when sim > theta1
				for(Clusters::iterator cluster_iter = pClusters->begin(); cluster_iter != pClusters->end(); cluster_iter++)
				{
					current_check++;
					if(current_check % 200 == 0)
						std::cout << "Adding/Skiping committee " << current_check << " / " << pClusters->size() << std::endl;
					Matrices::CompressedVectorPtr pCentroid = CountAvg( cluster_iter->second );
					if(!m_bUseTScoreRank)// for mi we hav got to weight the centroid... for rwf - avg is enough
						WeightCentroid(pCentroid);
					bool all_ok = true;
					// check all
					for(Committees::iterator ci = m_Committees.begin(); ci != m_Committees.end() && all_ok; ci++)
					{
						double sim = SimCentroid(pCentroid, (ci->first));
						if(sim > Config::theta1)
							all_ok = false;
					}
					//check added in this iteration
					for(Committees::iterator ci = current_committee.begin(); ci != current_committee.end() && all_ok; ci++)
					{
						double sim = SimCentroid(pCentroid, (ci->first));
						if(sim > Config::theta1)
							all_ok = false;
					}
					if(all_ok)
					{
						current_committee.push_back( std::pair<Matrices::CompressedVectorPtr, std::set<std::string> >
								(pCentroid, cluster_iter->second) );

						if (Config::verbose)
						{
							std::cout << "Adding committee: ";
							std::copy( cluster_iter->second.begin(), cluster_iter->second.end(),
									   std::ostream_iterator<std::string>(std::cout, " ") );
							std::cout << std::endl;
						}
					}
					else
					{
						if (Config::verbose)
						{
							std::cout << "Skipping committee: ";
							std::copy( cluster_iter->second.begin(), cluster_iter->second.end(),
									   std::ostream_iterator<std::string>(std::cout, " ") );
							std::cout << std::endl;
						}
					}
				}
				// step 4
				if(current_committee.size() != 0)
				{
					std::cout << "We have found " << current_committee.size() << " committees, searching for residues..." << std::endl;
					m_Committees.insert(m_Committees.end(), current_committee.begin(), current_committee.end() );
					std::vector<std::string> residues;

					//count residues (step5)
					for(size_t i = 0; i < words.size(); i++)
					{
						Matrices::CompressedVectorPtr v = m_pMatrix->GetRow( words[i] );
						bool all_lower = true;
						for(size_t j = 0; j < m_Committees.size() && all_lower; j++)
						{
							double sim = SimCentroid(v, m_Committees[j].first );
							if( sim > Config::theta2 )
								all_lower = false;
						}
						if(all_lower)
						{
							if(m_pMatrix->GetRowFeatures()->operator[](i).GetTermFrequency() > Config::committee_min_tf
									&& m_pMatrix->GetRowFeatures()->operator[](i).GetNZ() > Config::committee_min_nz )
							{
								if (Config::verbose)
									std::cout << "Adding to residues: " << words[i] << std::endl;

								residues.push_back(words[i]);
							}
						}

						if( i % 100 == 0)
							std::cout << "Searching for residues, word " << i << " / " << words.size() << std::endl;
					}
					words = residues;
				}
				else
				{
					std::cout << "No residues / current_committee is empty...\n";
					words.clear();
				}
				if (Config::intermediate)
					WriteCommittees(current_committee, iteration_number);

				if(words.size() )
					iteration_number ++;
			}//while we have something to do in phase 2 (recursive step)
		}//CBC::Phase2

		void CBC::WeightCentroid( Matrices::CompressedVectorPtr pCentroid)
		{
			bool discount = true;
			double row_tf = Measures::CountVectorSum(* pCentroid );
			for(std::map<size_t, double>::iterator iter = pCentroid->GetValues()->begin(); iter != pCentroid->GetValues()->end();iter++)
			{
				double val = iter->second;
				double column_tf = (double) (* m_pMatrix->GetColumnFeatures())[iter->first].GetTermFrequency();
				double denom = (row_tf + column_tf)/m_dN;
				double nom = val / m_dN;
				double mi = nom/denom;
				if( discount )
				{
					double a = val/(val+1);
					double m = row_tf <  column_tf? row_tf : column_tf;
					double b = m/(m+1);
					double discount =  a*b;
					mi *= discount;
					pCentroid->SetElement( iter->first, mi);
				}
				else pCentroid->SetElement( iter->first, mi);
			}
		}

		double CBC::SimCentroid(Matrices::CompressedVectorPtr &v1, Matrices::CompressedVectorPtr &v2)
		{
			static smartcomparator::simfunctions::LocalCosineSimilarity sim(false, false);

			return sim.CountSimilarity(*v1,*v2);
		}

		void CBC::Phase3()
		{
			if (Config::remove_features == WEIGHTED)
				CountFeaturesSums();

			using namespace boost::posix_time;
			std::cout << "Phase 3... Assigning words to committees\n";
			// for each word e
			ptime start = microsec_clock::local_time();
			time_duration total_time;
			for(size_t word_no = 0; word_no < m_pMatrix->NumberOfRows(); word_no++ )
			{
				// (count similarities for each word)
				// S <- top-200 most similar committees to e, sorted by similarity to e
				std::multimap<double, std::pair<Committee*, size_t>, std::greater<double> > S; // sim -> (Committee, id in m_Committees)
				const std::string& word = m_pMatrix->GetRowFeatures()->operator[](word_no).GetLabel();
				Matrices::CompressedVectorPtr word_vector_copy(new Matrices::CompressedVector(*( m_pMatrix->GetRow(word)) ));
				PrepareCommitteesForWord( S, word);
				// C = output clusters of e
				std::vector<Committee*> C;
				// while !empty( S ):
				while( S.size() )
				{
					// c = S.pop() ; c - most similar commitee
					std::multimap<double, std::pair<Committee*, size_t>, std::greater<double> >::iterator rbegin = S.begin();
					size_t no = rbegin->second.second;
					if( no > m_Committees.size() )
						std::cerr << "1no > committees size!!! csize=" << C.size() << " s_size=" << S.size() << std::endl;
					Committee* c = rbegin->second.first;
					Matrices::CompressedVectorPtr c_vector = c->first;
					no = rbegin->second.second;
					double sim = SimCentroid( c_vector, word_vector_copy);
					if( sim < Config::sigma )
						break;
					if ( C.size() == 0 || !IsSimilar( c, C) )
					{
						size_t no = rbegin->second.second;
						WordSimPair wsp;
						wsp.word = word;
						wsp.sim = sim;
						m_mWordToCommitteeIndex[word].insert( no );
						m_mCommitteeIndexToWord[ no ].insert(wsp);
						RemoveOverlapingFeatures(word_vector_copy, c->first);
						C.push_back(c);
						S.erase(rbegin);
						RecalculateSimilariteis( S, word);
					}
					else
						S.erase(rbegin);

				}//while S.size() (a.k.a. !empty(S) )

				ptime end = microsec_clock::local_time();
				time_duration len = end - start;
				if( len.total_milliseconds() > 5000)
				{
					total_time += len;
					long total_time_ms = total_time.total_milliseconds();
					double words_per_sec = double( (word_no+1)*1000)  / total_time_ms;
					long eta_ms = (m_pMatrix->NumberOfRows()*total_time_ms) / (word_no+1) - total_time_ms;
					std::cout << "Assigning words to committees, " << word_no << " / " << m_pMatrix->NumberOfRows()
						<< "| " << " words/sec=" << words_per_sec << " eta: " << to_simple_string(milliseconds(eta_ms) )
						<< std::endl;
					start = microsec_clock::local_time();
				}
			}
		}//CBC::Phase3

		bool CBC::IsSimilar(Committee* c, std::vector<Committee*> &all)
		{
			Matrices::CompressedVectorPtr c_v = c->first;
			for(size_t i =  0; i < all.size(); i++)
			{
				double sim = SimCentroid(c_v, all[i]->first);
				if( sim < Config::not_similiar_threshold )
					return false;
			}
			return true;
		}//bool IsNotSimilar(Committee* c, std::vector<Committee*> &all);

		void CBC::RemoveOverlapingFeatures(Matrices::CompressedVectorPtr word_vector_copy, Matrices::CompressedVectorPtr c)
		{
			//FIXME: jako, ze uzywamy kosinusa to mozemy tylko iterowac po indeksach w centroidzie... ale w ogolnosci to jest niepoprawne!

			std::map<size_t, double> *word_vals = word_vector_copy->GetValues();
			std::map<size_t, double> *centroid_vals = c->GetValues();


			if (Config::remove_features == NORMAL)
			{
				for(std::map<size_t, double>::iterator iter = centroid_vals->begin(); iter != centroid_vals->end(); iter++)
					word_vals->erase( iter->first );
			}
			if (Config::remove_features == MIN)
			{
				for(std::map<size_t, double>::iterator iter = centroid_vals->begin(); iter != centroid_vals->end(); iter++)
				{
					double val_committee = iter->second;
					std::map<size_t, double>::iterator find = word_vals->find( iter->first );
					if ( find == word_vals->end() )
						word_vals->erase(iter->first);
					else
					{
						double val_vector = find->second;
						if( val_committee < val_vector)
							find->second = val_committee;
					}
				}
			}
			if (Config::remove_features == WEIGHTED)
			{
				for(std::map<size_t, double>::iterator iter = centroid_vals->begin(); iter != centroid_vals->end(); iter++)
				{
					double val_committee = iter->second;
					std::map<size_t, double>::iterator find = word_vals->find( iter->first );
					if ( find == word_vals->end() )
						word_vals->erase(iter->first);
					else
					{
						double val_vector = find->second;
						double w = val_vector - val_vector* (val_committee/m_vSumCols[iter->first] );
						find->second = w;

					}
				}
			}
		}

		void CBC::CountFeaturesSums()
		{
			std::cout << "Counting features sums for committees..." << std::endl;
			for(size_t col_no = 0; col_no < m_pMatrix->NumberOfColumns(); col_no++)
			{
				double sum = 0.0;
				for(Committees::iterator iter = m_Committees.begin(); iter != m_Committees.end(); iter++)
				{
					Committee &c = *iter;
					Matrices::CompressedVectorPtr pCenteroid = c.first;
					sum += (*pCenteroid)[ col_no ];
				}
				m_vSumCols.push_back( sum );

				if(col_no % 50 )
					std::cout << "Counting features sums for committees..." << col_no << " / " << m_pMatrix->NumberOfColumns() << std::endl;

			}
			std::cout << "Counting features sums...Done" << std::endl;
		}

		void CBC::RecalculateSimilariteis(std::multimap<double,std::pair<Committee*, size_t>, std::greater<double> > &sims, const std::string &word)
		{
			std::vector< std::pair<Committee*, size_t> > all;
			for(std::multimap<double, std::pair<Committee*, size_t> >::iterator iter = sims.begin(); iter != sims.end() ; iter++)
				all.push_back( iter->second );
			sims.clear();
			Matrices::CompressedVectorPtr v = m_pMatrix->GetRow(word);
			for(size_t i = 0; i < all.size(); i++)
			{
				std::pair<Committee *, size_t> &c = all[i];
				double sim = SimCentroid( v, c.first->first );
				sims.insert(std::pair< double, std::pair<Committee*, size_t> > (sim, c) );
			}
		}

		void CBC::PrepareCommitteesForWord(std::multimap<double,std::pair<Committee*,size_t>, std::greater<double> > &out, const std::string &word)
		{
			out.clear();
			Matrices::CompressedVectorPtr v = m_pMatrix->GetRow(word);
			// count sims for all
			for(size_t committee_no = 0; committee_no < m_Committees.size(); committee_no++)
			{
				double sim = SimCentroid( v, m_Committees[committee_no].first );
				Committee * c = & m_Committees[committee_no];
				out.insert( std::pair<double, std::pair< Committee*, size_t> >( sim, std::pair<Committee*, size_t>(c, committee_no ) ) );
			}
			// use only first X
			size_t added = 0;
			for(std::multimap<double, std::pair<Committee*, size_t> >::iterator iter = out.begin(); iter != out.end() ; iter++)
			{
				if(added++ >= Config::phase_3_top_x)
				{
					out.erase(iter, out.end() );
					break;
				}
			}
		}//PrepareCommitteesForWord

		Matrices::CompressedVectorPtr CBC::CountAvg( std::set<std::string> &cluster)
		{
			Matrices::CompressedVectorPtr v( new Matrices::CompressedVector(m_pMatrix->NumberOfColumns() ) );
			size_t n = 0;
			for( std::set<std::string>::iterator iter = cluster.begin(); iter != cluster.end(); iter++)
			{
				Matrices::CompressedVectorPtr part = m_pMatrix->GetRow( *iter );
				if(part == NULL)
				{
					std::cerr << "Warning! No row in matrix " << m_pMatrix->GetName() << " for word: " << *iter << ". Skipping...\n";
					continue;
				}
				n++;
				std::map<size_t, double> * vals = part->GetValues();
				for(std::map<size_t, double>::iterator i = vals->begin(); i != vals->end(); i++)
				{

					double val = i->second + (*v)[i->first];
					if(val > 0.00000001 || val <-0.00000001)
						v->SetElement(i->first, val);
				}
			}
			std::map<size_t, double> * vals = v->GetValues();
			std::set<size_t> toerase;
			for(std::map<size_t, double>::iterator i = vals->begin(); i != vals->end(); i++)
			{
				double val = i->second / (double)n;
				if(val > 0.00000001 || val < -0.00000001)
					i->second = val;
				else
					toerase.insert(i->first);
			}
			for(std::set<size_t>::iterator i = toerase.begin(); i != toerase.end(); i++)
				vals->erase(*i);
			return v;

		}//CountAvg

		void CBC::WriteCommittees(Committees &committees, size_t iter_no)
		{
			boost::filesystem::path p= Config::workdir;
			p /= "committees." + boost::lexical_cast<std::string>(iter_no);

			if (Config::verbose)
				std::cout << "Writing committees to: " << p.string() << std::endl;

			std::ofstream file( p.string().c_str() );
			for(Committees::iterator ci = committees.begin(); ci != committees.end(); ci++)
			{
				std::copy( ci->second.begin(), ci->second.end(), std::ostream_iterator<std::string>(file, ";") );
				file << ": ";
				std::map<size_t, double> *vals = ci->first->GetValues();
				for(std::map<size_t, double>::iterator iter = vals->begin(); iter != vals->end(); iter++)
					file << iter->first << " " << iter->second << " ";
				file << std::endl;
			}
		}

		void CBC::ResumePhase3(std::string & prev_phase_filename)
		{
			std::cout << "Reading committees from file: " << prev_phase_filename << std::endl;
			std::ifstream file(prev_phase_filename.c_str() );
			if( !file )
				throw std::runtime_error( "Cannot open file: " + prev_phase_filename );
			std::string line;
			//int line_no  = 0;
			while( std::getline(file,line) )
			{
				boost::trim(line);
				std::vector<std::string> first_split;
				boost::split(first_split, line, boost::is_any_of(":") );
				if( first_split.size() != 2 )
				{
					std::cerr << "Wrong line: " << line << std::endl;
					continue;
				}
				std::vector<std::string> vWords;
				boost::split(vWords, first_split[0], boost::is_any_of(";") );
				std::set<std::string> words;
				for(size_t i = 0; i < vWords.size(); i++)
				{
					std::string w = vWords[i];
					boost::trim(w);
					if(w.length() )
						words.insert(w);
				}
				if(words.size() == 0 )
				{
					std::cerr << "Wrong line: " << line << std::endl;
					continue;
				}
				Matrices::CompressedVectorPtr v(new Matrices::CompressedVector( m_pMatrix->NumberOfColumns() ) );
				std::vector<std::string> strVals;
				std::string line2 = first_split[1];
				boost::trim(line2);
				boost::split(strVals, line2, boost::is_any_of(" ") );
				bool ok = false;

				for(size_t i = 0; i < strVals.size() - 1; i+=2 )
				{
					ok=true;
					size_t id = boost::lexical_cast<size_t>( strVals[i] );
					double val = boost::lexical_cast<double>(strVals[i+1] );
					v->SetElement(id, val);
				}
				if(ok)
				{
					m_Committees.push_back( Committee( v, words ) );
				}

			}
			std::cout << "Read " << m_Committees.size() << " committees..." << std::endl;
			Phase3();
			WriteResults();

		}//void ResumePhase3(std::string & prev_phase_filename);
	}
}//clustering

