/*
 * AresV.cpp
 *
 *  Created on: 11-Jul-08
 *      Author: indurks
 */

#include "AresV.h"
#include <iomanip>

//#define BASIC_INPUT_DEBUG

void HandleChemical(int c, GuardedTables &myGuardedTables, ChemicalMatrix &myCM, list<MyCloud> &myClouds) {
	bool verbose = false;

	if(verbose)
		cout << "Creating a cloud for chemical: " << c << std::endl;

	MyCloud myCloud(c);
	myClouds.push_back(myCloud);

	// Go through and make a list of chemicals that cross the sub tree threshold
	list<int> subTreeCandidates;

	for(map<int, list<int> >::iterator iter = myCM.chemMatrix[c].second.begin(); iter != myCM.chemMatrix[c].second.end(); iter++)
		if(iter->second.size() >= SUB_TREE_THRESHOLD)
		{
			if(!(iter->first == 0 && NULL_SUBTREES_ALLOWED == false))
			{
				if(verbose)
					cout << "\tCreating a subTree for chemical " << iter->first << " with " << iter->second.size() << " entries" << std::endl;

				subTreeCandidates.push_back(iter->first);
			}
		}

	// Then make the sub-trees, and then add them to the cloud
	for(list<int>::iterator iter = subTreeCandidates.begin(); iter != subTreeCandidates.end(); iter++)
	{
		vector<Reaction> rxns;

		for(list<int>::iterator jiter = myCM.chemMatrix[c].second[*iter].begin(); jiter != myCM.chemMatrix[c].second[*iter].end(); jiter++)
			rxns.push_back(myGuardedTables.getRxn(*jiter));

		for(vector<Reaction>::iterator jiter = rxns.begin(); jiter != rxns.end(); jiter++)
			myCM.RemoveReaction(*jiter);

		myClouds.back().AddTree(new MyTree(*iter, rxns, *iter));
	}

	if(verbose)
		cout << "\tNow Forming Main Tree!" << std::endl;

	// Then make the main tree, and add it to the cloud.
	vector<Reaction> mainTreeRxns;

	for(map<int, list<int> >::iterator iter = myCM.chemMatrix[c].second.begin(); iter != myCM.chemMatrix[c].second.end(); iter++)
		for(list<int>::iterator jiter = iter->second.begin(); jiter != iter->second.end(); jiter++)
			mainTreeRxns.push_back(myGuardedTables.getRxn(*jiter));

	for(vector<Reaction>::iterator iter = mainTreeRxns.begin(); iter != mainTreeRxns.end(); iter++)
		myCM.RemoveReaction(*iter);

	myClouds.back().SetMainTree(new MyTree(0, mainTreeRxns, c));
	myClouds.back().CalculateSumPropensity();
	myClouds.back().SortCloud();

	if(myClouds.back().mainTree->myRxns.size() == 0 && myClouds.back().trees.size() ==0)
		myClouds.pop_back();
}

enum VM_INSTRUCTION {
	UPDATE_SINGLE_RXN_IN_SUPERCACHE,	// 0
	UPDATE_DOUBLET_RXN_IN_SUPERCACHE,	// 1
	SET_CLOUD,							// 2
	UPDATE_RXN_IN_CLOUD,				// 3
	UPDATE_SUBTREE_IN_CLOUD,			// 4
	UPDATE_CLOUD_NO_PRIMARY,			// 5
	UPDATE_CLOUD_WITH_PRIMARY,			// 6
	SET_AND_UPDATE_CLOUD_ONE_PASS,		// 7
	SC_URIC_UCNP,						// 8
	SC_URIC_UCNP_UP_ONE,				// 9
	SC_URIC_UCNP_UP_TWO,				// 10
	SC_URIC_UCNP_UP_THREE				// 11
};

class ChemicalInstruction {
	int chem;
	vector<bool> touchedClouds;
	map<int, list< pair<int, int> > > instructions;
	list<pair<int, int> > superCacheInstructions;

public:
	list< pair<int, int> > finalInstructions;
	list< list<int> > realFinalInstructions;

	ChemicalInstruction(int chemical, int numClouds) {
		chem = chemical;

		for(int i = 0; i < numClouds; i++)
			touchedClouds.push_back(false);
	}

	int getChem() {
		return chem;
	}

	void AddInstructionForSuperCache(pair<int, int> p) {
		superCacheInstructions.push_back(p);
	}

	void AddInstructionForCloud(int cloud, pair<int, int> p) {
		TagCloudForOpening(cloud);
		instructions[cloud].push_back(p);
	}

	void TagCloudForOpening(int cloud) {
		touchedClouds[cloud] = true;
	}

	void WrapUpChemInstruction(vector<int> &myCloudChemMapping) {
		cout << "Wrapping up Chem Instruction for chem = " << chem << " and touchedClouds.size() = " << touchedClouds.size() << std::endl;

		// Go through the super-cached reactions
		for(list<pair<int, int> >::iterator iter = superCacheInstructions.begin(); iter != superCacheInstructions.end(); iter++)
			finalInstructions.push_back(*iter);

		cout << "\tFinished Wrapping super-cached reactions..." << std::endl;

		// Go through each of the touched trees
		for(int i = 0; i < (signed)touchedClouds.size(); i++) {
			if(touchedClouds[i]) {
				cout << "\tGot through for i = " << i << std::endl;
				cout << "myCloudChemMapping.size() = " << myCloudChemMapping.size() << std::endl;
				cout << "myCloudChemMapping[i] = " << myCloudChemMapping[i] << std::endl;

				if(chem == myCloudChemMapping[i]) {
					// See if it's going to be a simple SET_AND_UPDATE_TREE_ONE_PASS
					if(instructions[i].size() == 0) {
						// Sure is!
						pair<int, int> p(SET_AND_UPDATE_CLOUD_ONE_PASS, i);
						finalInstructions.push_back(p);
					} else {
						instructions[i].sort();
						pair<int, int> p1(SET_CLOUD, i);
						pair<int, int> p2(UPDATE_CLOUD_WITH_PRIMARY, i);

						finalInstructions.push_back(p1);
						for(list< pair<int, int> >::iterator iter = instructions[i].begin(); iter != instructions[i].end(); iter++)
							finalInstructions.push_back(*iter);
						finalInstructions.push_back(p2);
					}
				} else {
					if(instructions[i].size() > 0)
					{
						instructions[i].sort();
						pair<int, int> p1(SET_CLOUD, i);
						pair<int, int> p2(UPDATE_CLOUD_NO_PRIMARY, i);

						finalInstructions.push_back(p1);
						for(list< pair<int, int> >::iterator iter = instructions[i].begin(); iter != instructions[i].end(); iter++)
							finalInstructions.push_back(*iter);
						finalInstructions.push_back(p2);
					}
					else
					{
						cout << "LOLOLOL" << std::endl;
						exit(0);
					}
				}
			}
		}

		// Now we want to detect and catch 2-3-5's
		vector<pair<int, int> > GFIB;

		for(list< pair<int, int> >::iterator iter = finalInstructions.begin(); iter != finalInstructions.end(); iter++)
			GFIB.push_back(*iter);

		finalInstructions.clear();

		while(!GFIB.empty())
		{
			list<int> zebracat;
			if(GFIB.size() >= 3) {
				if(GFIB[0].first == SET_CLOUD && GFIB[1].first == UPDATE_RXN_IN_CLOUD && GFIB[2].first == UPDATE_CLOUD_NO_PRIMARY)
				{
					// Yay now we can do an optimization!!!
					zebracat.push_back(SC_URIC_UCNP);
					zebracat.push_back(GFIB[0].second);
					zebracat.push_back(GFIB[1].second);
					GFIB.erase(GFIB.begin());
					GFIB.erase(GFIB.begin());
					GFIB.erase(GFIB.begin());
				} else {
					zebracat.push_back(GFIB[0].first);
					zebracat.push_back(GFIB[0].second);
					GFIB.erase(GFIB.begin());
				}
			} else {
				zebracat.push_back(GFIB[0].first);
				zebracat.push_back(GFIB[0].second);
				GFIB.erase(GFIB.begin());
			}

			realFinalInstructions.push_back(zebracat);
		}
	}
};

/* PARAMETERS ARE:
	1) Main Rxn input file
	2) Main Chem input file
	3) Sum-Prop file
	4) Super-Cache Cutoff file
	5) Output Header
	6) MTDF_OUTPUT*/
int main(int argc, char* argv[]) {
	cout << setprecision(10);

	vector<Reaction> myRxns;
	vector<Chemical> myChems;
	list<MyCloud> myClouds;

	// Load the chem matrix in
	cout << "Loading Chemical Matrix and Guarded Table...	";

	ImportChemsAndReactionsAndSumProp(string(argv[1]), string(argv[2]), string(argv[3]), myRxns, myChems);

	GuardedTables myGuardedTables(myRxns, myChems);

	cout << "COMPLETE" << std::endl;
	cout << "Read in " << myRxns.size() << " reactions." << std::endl;
	cout << "Read in " << myChems.size() << " species." << std::endl;

	ChemicalMatrix myCM(myGuardedTables);

	// Create our super-cache
	cout << std::endl << "Forming Super-Cache...	";

	ifstream SC_Cutoff_Input(argv[4]);
	int numSuperCachedReactions;
	SC_Cutoff_Input >> numSuperCachedReactions;

	vector<Reaction> SuperCachedRxns;

	for(int i = 0; i < numSuperCachedReactions; i++)
		SuperCachedRxns.push_back(myCM.RemoveReaction(myGuardedTables.getRxn(i)));

	cout << "COMPLETE" << std::endl;
	cout << "Super-Cached " << numSuperCachedReactions << " reactions." << std::endl << std::endl;

	// Create a queue of chemicals to go through
	cout << "Creating a queue of chemicals to process...	";

	list<pair<double, int> > queueOfChems;

	for(SuperDuperMap::iterator iter = myCM.chemMatrix.begin(); iter != myCM.chemMatrix.end(); iter++) {
		if(iter->first != 0)
		{
			pair<double, int> p(myCM.getScore(iter->first), iter->first);
			queueOfChems.push_back(p);
		}
	}

	queueOfChems.sort();

	if(myCM.hasChemical(0))
	{
		pair<double, int> p(myCM.getScore(0), 0);
		queueOfChems.push_front(p);
	}

	cout << "COMPLETE" << std::endl;
	cout << "Number of Chemicals in Cheme Queue = " << queueOfChems.size() << std::endl << std::endl;

	// For each chemical in the queue, pass it to Handle Chemical
	cout << "Handling chemicals in queue...	";

	/*for(list<pair<double, int> >::iterator iter = queueOfChems.begin(); iter != queueOfChems.end(); iter++)
		cout << "\tChem " << iter->second << " has score " << iter->first << std::endl;*/

	for(list<pair<double, int> >::reverse_iterator iter = queueOfChems.rbegin(); iter != queueOfChems.rend(); iter++) {
		if(myCM.chemMatrix.find(iter->second) != myCM.chemMatrix.end())
		{
			HandleChemical(iter->second, myGuardedTables, myCM, myClouds);
		}
	}

	cout << "COMPLETE" << std::endl;
	cout << "Number of Clouds = " << myClouds.size() << std::endl << std::endl;

	// Build the main tree
	cout << "Constructing Main Tree...	" << std::endl;
	int myB, myA, myD;
	InitializeTree(myB, myA, myD, myClouds.size());

	int MainTreeOffset = myB;
	int SizeOfMainTree = myA;
	int DepthOfMainTree = myD;
	cout << "COMPLETE" << std::endl << std::endl;

	// Generate the instruction table
	cout << "Generating Instruction Table...	" << std::endl;

	vector<ChemicalInstruction> InstructionTable;

	for(int i = 0; i < (signed)myChems.size(); i++)
	{
		ChemicalInstruction ci(i, myClouds.size());
		InstructionTable.push_back(ci);
	}

	// Go through the Super-Cache
	cout << "\tProcessing Super-Cache...	";

	for(int i = 0; i < (signed)SuperCachedRxns.size(); i++)
		if(SuperCachedRxns[i].r1 == SuperCachedRxns[i].r2) {
			// We have a doublet
			pair<int, int> p(UPDATE_DOUBLET_RXN_IN_SUPERCACHE, i);
			if(SuperCachedRxns[i].r1 != 0)
				InstructionTable[SuperCachedRxns[i].r1].AddInstructionForSuperCache(p);
		} else {
			pair<int, int> p(UPDATE_SINGLE_RXN_IN_SUPERCACHE, i);
			if(SuperCachedRxns[i].r1 != 0)
				InstructionTable[SuperCachedRxns[i].r1].AddInstructionForSuperCache(p);
			if(SuperCachedRxns[i].r2 != 0)
				InstructionTable[SuperCachedRxns[i].r2].AddInstructionForSuperCache(p);
		}

	cout << "COMPLETE HELLO" << std::endl;
	cout << "\tProcessing Clouds...	" << std::endl;

	int zebra = 0;

	for(list<MyCloud>::iterator iter = myClouds.begin(); iter != myClouds.end(); iter++, zebra++)
	{
		cout << "\t\t" << "Processing Cloud with chemical = " << iter->deMultiplierChemical << std::endl;

		if(iter->deMultiplierChemical != 0)
			InstructionTable[iter->deMultiplierChemical].TagCloudForOpening(zebra);

		// First go through the main tree
		for(int i = 0; i < (signed)iter->mainTree->myRxns.size(); i++)
		{
			pair<int, int> p(UPDATE_RXN_IN_CLOUD, i);
			int affectedChem = iter->mainTree->myRxns[i].r2;
			if(affectedChem != 0)
				InstructionTable[affectedChem].AddInstructionForCloud(zebra, p);
		}

		// Now go through each subtree in the cloud
		for(int i = 0; i < (signed)iter->trees.size(); i++)
		{
			pair<int, int> p(UPDATE_SUBTREE_IN_CLOUD, i);
			int affectedChem = iter->trees[i]->deMultiplierChemical;
			if(affectedChem != 0)
				InstructionTable[affectedChem].AddInstructionForCloud(zebra, p);
		}
	}

	vector<int> myCloudChemMapping;

	for(list<MyCloud>::iterator iter = myClouds.begin(); iter != myClouds.end(); iter++)
		myCloudChemMapping.push_back(iter->deMultiplierChemical);

	cout << "About to wrap up instructions!" << std::endl;

	for(vector<ChemicalInstruction>::iterator iter = InstructionTable.begin(); iter != InstructionTable.end(); iter++)
		iter->WrapUpChemInstruction(myCloudChemMapping);

	cout << "COMPLETE" << std::endl << std::endl;

	cout << "Generating Output Instruction Table... " << std::endl;

	string outputHeader(argv[5]);
	string strBasicOutputFileName = outputHeader + "_BASIC_OUTPUT.dat";
	string strRxnOutputFileName = outputHeader + "_RXN_OUTPUT.dat";
	string strSpeciesOutputFileName = outputHeader + "_SPECIES_OUTPUT.dat";

	ofstream basicOutput(strBasicOutputFileName.c_str());
	ofstream rxnOutput(strRxnOutputFileName.c_str());
	ofstream speciesOutput(strSpeciesOutputFileName.c_str());

#ifdef BASIC_INPUT_DEBUG
	basicOutput << "$%	numChemicals	%$" << std::endl;
#endif
	basicOutput << myGuardedTables.getNumChems() << std::endl;

#ifdef BASIC_INPUT_DEBUG
	basicOutput << "$%	numReactions	%$" << std::endl;
#endif
	basicOutput << myGuardedTables.getNumRxns() << std::endl;

#ifdef BASIC_INPUT_DEBUG
	basicOutput << "$%	num-Super-Cached-Reactions	%$" << std::endl;
#endif
	basicOutput << numSuperCachedReactions << std::endl;

#ifdef BASIC_INPUT_DEBUG
	basicOutput << "$%	numClouds	%$" << std::endl;
#endif
	basicOutput << myClouds.size() << std::endl;

#ifdef  BASIC_INPUT_DEBUG
	basicOutput << "$%	MainTreeOffset	%$" << std::endl;
#endif
	basicOutput << MainTreeOffset << std::endl;

	int numRxnsOutputted = 0;

	vector<int> rxnTape, rxnTapeIndex;

	cout << "\tOutputting Super Cached Reactions... ";

	// Now output species
	for(int i = 0; i < myGuardedTables.getNumChems(); i++)
		speciesOutput << myGuardedTables.getChem(i).initConc << std::endl;

	// Now output super cached reactions
	for(int i = 0; i < SuperCachedRxns.size(); i++)
	{
		list<int> stuff = getRxnTapeFromReaction(SuperCachedRxns.at(i));
		rxnOutput << SuperCachedRxns[i].r1 << "	" << SuperCachedRxns[i].r2 << "	" << SuperCachedRxns[i].p1 << "	" << SuperCachedRxns[i].p2 << "	" << SuperCachedRxns[i].rate << std::endl;

		rxnTape.push_back(stuff.size()/2);

		for(list<int>::iterator iter = stuff.begin(); iter != stuff.end(); iter++)
			rxnTape.push_back(*iter);
	}

	numRxnsOutputted += SuperCachedRxns.size();

	cout << "COMPLETE." << std::endl;

	// Now comes the tricky part where we have to simultaneously emit both reactions and instructions
	for(list<MyCloud>::iterator iter = myClouds.begin(); iter != myClouds.end(); iter++)
	{
		cout << "\tOutputting cloud for chemical " << iter->deMultiplierChemical << "... ";
#ifdef BASIC_INPUT_DEBUG
		basicOutput << "$%	Cloud for Chemical " << iter->deMultiplierChemical << "	%$" << std::endl;
		basicOutput << "$%	-------------------------------	%$" << std::endl;
		basicOutput << "$%	deMultiplierChemical	%$" << std::endl;
#endif
		basicOutput << iter->deMultiplierChemical << std::endl;

#ifdef BASIC_INPUT_DEBUG
		basicOutput << "$%	mainTreeOffset	%$" << std::endl;
#endif
		basicOutput << iter->mainTree->myBefore << std::endl;

#ifdef BASIC_INPUT_DEBUG
		basicOutput << "$%	mainTreeReboost	%$" << std::endl;
#endif
		basicOutput << numRxnsOutputted << std::endl;

#ifdef BASIC_INPUT_DEBUG
		basicOutput << "$%	numSubTrees	%$" << std::endl;
#endif
		basicOutput << iter->trees.size() << std::endl;

#ifdef BASIC_INPUT_DEBUG
		basicOutput << "$%	mainTree-NumReactions	%$" << std::endl;
#endif
		basicOutput << iter->mainTree->numRxns << std::endl;

		// Now output the reactions
		for(int i = 0; i < iter->mainTree->myRxns.size(); i++) {
			rxnOutput << iter->mainTree->myRxns.at(i).r1 << "	" << iter->mainTree->myRxns.at(i).r2 << "	" << iter->mainTree->myRxns.at(i).p1 << "	" << iter->mainTree->myRxns.at(i).p2 << "	" << iter->mainTree->myRxns.at(i).rate << std::endl;

			list<int> stuff = getRxnTapeFromReaction(iter->mainTree->myRxns.at(i));
			rxnTape.push_back(stuff.size() / 2);

			for(list<int>::iterator jiter = stuff.begin(); jiter != stuff.end(); jiter++)
				rxnTape.push_back(*jiter);
		}

		numRxnsOutputted += iter->mainTree->myRxns.size();

		// Now output the subTrees
		for(int i = 0; i < iter->trees.size(); i++) {
			MyTree *tempTree = iter->trees[i];

			basicOutput << std::endl;

#ifdef BASIC_INPUT_DEBUG
			basicOutput << "\t" << "$% subTree for 2nd Chemical: " << tempTree->deMultiplierChemical << "	%$" << std::endl;
			basicOutput << "\t" << "$% --------------------------------- %$" << std::endl;
			basicOutput << "\t" << "$% offset %$" << std::endl;
#endif
			basicOutput << "\t" << tempTree->myBefore << std::endl;

#ifdef BASIC_INPUT_DEBUG
			basicOutput << "\t" << "$% reBoost %$" << std::endl;
#endif
			basicOutput << "\t" << numRxnsOutputted << std::endl;

#ifdef BASIC_INPUT_DEBUG
			basicOutput << "\t" << "$% deMultiplierChemical %$" << std::endl;
#endif
			basicOutput << "\t" << tempTree->deMultiplierChemical << std::endl;

#ifdef BASIC_INPUT_DEBUG
			basicOutput << "\t" << "$% numRxns %$" << std::endl;
#endif
			basicOutput << "\t" << tempTree->numRxns << std::endl;
			basicOutput << std::endl;

			// Now output the reactions
			for(int j = 0; j < tempTree->numRxns; j++)
			{
				rxnOutput << tempTree->myRxns[j].r1 << "	" << tempTree->myRxns[j].r2 << "	" << tempTree->myRxns[j].p1 << "	" << tempTree->myRxns[j].p2 << "	" << tempTree->myRxns[j].rate << std::endl;

				list<int> stuff = getRxnTapeFromReaction(tempTree->myRxns[j]);

				rxnTape.push_back(stuff.size() / 2);

				for(list<int>::iterator jiter = stuff.begin(); jiter != stuff.end(); jiter++)
					rxnTape.push_back(*jiter);
			}

			numRxnsOutputted += tempTree->numRxns;
		}

		basicOutput << std::endl;

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

	basicOutput << std::endl << std::endl;

	for(vector<ChemicalInstruction>::iterator iter = InstructionTable.begin(); iter != InstructionTable.end(); iter++)
	{
#ifdef BASIC_INPUT_DEBUG
		basicOutput << "$% Num Instructions for chemical: " << iter->getChem() << " %$" << std::endl;
#endif
		basicOutput << iter->realFinalInstructions.size() << std::endl;

#ifdef BASIC_INPUT_DEBUG
		basicOutput << "\t$% --- %$" << std::endl;
#endif
		for(list<list<int> >::iterator jiter = iter->realFinalInstructions.begin(); jiter != iter->realFinalInstructions.end(); jiter++)
		{
			basicOutput << "\t";
			for(list<int>::iterator kiter = jiter->begin(); kiter != jiter->end(); kiter++)
				basicOutput << (*kiter) << "	";
			basicOutput << std::endl;
		}
	}

	basicOutput << std::endl;

	for(int i = 0, k = 0; k < myGuardedTables.getNumRxns(); k++)
	{
#ifdef BASIC_INPUT_DEBUG
		basicOutput << "$% Num Entries for Reaction: " << k << " %$" << std::endl;
#endif
		basicOutput << rxnTape.at(i) << std::endl;
		int numEntriesOnTape = rxnTape.at(i);
		i += 1;
		for(int j = 0; j < numEntriesOnTape; j++, i += 2)
			basicOutput << "\t" << rxnTape.at(i) << "	" << rxnTape.at(i + 1) << std::endl;
	}

	// Now do a basic check. The number of rxnsOutputted should be the same as the number of reactions
	if(numRxnsOutputted == myGuardedTables.getNumRxns())
		cout << "Successfully outputted all reactions in system." << std::endl;

	cout << "myClouds.size() = " << myClouds.size() << std::endl;
	cout << "Outputting MTDF... = ";
	
	if(1 == (DepthOfMainTree - 1)) {
		ofstream outputMTDF("ONE_MTDF.dat");
		outputMTDF << "ONE" << std::endl;
		cout << "ONE" << std::endl;
	}
	else if(2 == (DepthOfMainTree - 1)) {
		ofstream outputMTDF("TWO_MTDF.dat");
		outputMTDF << "TWO";
		cout << "TWO" << std::endl;
	}
	else if(3 == (DepthOfMainTree - 1)) {
		ofstream outputMTDF("THREE_MTDF.dat");
		outputMTDF << "THREE";
		cout << "THREE" << std::endl;
	}
	else if(4 == (DepthOfMainTree - 1)) {
		ofstream outputMTDF("FOUR_MTDF.dat");
		outputMTDF << "FOUR";
		cout << "FOUR" << std::endl;
	}
	else {
		ofstream outputMTDF("PLAIN_MTDF.dat");
		outputMTDF << "PLAIN";
		cout << "PLAIN" << std::endl;
	}

	cout << "PROGRAM COMPLETE." << std::endl;

	return 0;
}
