#ifndef _AmpToolsConfigGenerator_
#define _AmpToolsConfigGenerator_

#include <deque>
#include <map>
#include <set>
#include <vector>
#include <string>
#include <vector>
#include <sstream>
#include <fstream>
#include <iostream>
#include <stdlib.h>
#include <unistd.h>

#include "particleType.h"
#include "HCA_Amplitudes/resonanceinfo.h"

#include "HCA_TensorContraction/HCA_TensorContractor.h"

#include "Particle.h"
#include "ResonanceDecay.h"
#include "Transition.h"
#include "FitParamInfo.h"

using namespace std;

class HCA_ConfigGenerator
{
	public:
		HCA_ConfigGenerator(Particle_t locBeamPID, Particle_t locTargetPID, deque<Particle_t> locFinalStatePIDs) : 
		dBeamPID(locBeamPID), dTargetPID(locTargetPID), dFinalStatePIDs(locFinalStatePIDs)
		{
			dAmplitudeCounter = 0;
			dBinInSqrtSFlag = false;
			dDebugFlag = false;
			dNoBreitWignersFlag = false;
			dMaxNumFitIterations = 5000;

			dLimitToOneContractionSchemeFlag = true;
		}

		//Beam Energy //For MC Generation
		double dBeamEnergy; //In Lab Frame

		//PIDs
		Particle_t dBeamPID;
		Particle_t dTargetPID;
		deque<Particle_t> dFinalStatePIDs;

		//Breit-Wigner Control Params
		bool dBinInSqrtSFlag;
		bool dNoBreitWignersFlag;

		//Fit Control Params
		int dMaxNumFitIterations;

		//Limit Flags
		bool dLimitToOneContractionSchemeFlag;

		string dFitName;

		//data file names
		string dGenerateDataFileName;
		string dAcceptedDataFileName;
		string dExperimentDataFileName;

		//other file names
		string dConfigFileName;
		string dSeedFileName;
		string dNormalizationIntegralFileName;

		bool dDebugFlag;
		bool dIncludePhaseSpaceFlag;

		void Set_Particles(Particle* locPsiBeam, Particle* locPsiTarget, deque<Particle*>& locFinalStateParticles);
		void Set_FitParams(set<FitParamInfo*>& locFitParams);

		void Generate(deque<NodeBase*>& locTransitionTrees);
		void Print_Tree(NodeBase* locTreeBase, int locNumTabs = 0);
		void Print_ResonanceNode(ResonanceNode* locCurrentNode, int locNumTabs = 0);

		deque<ResonanceDecay*> dResonanceDecays;

	private:

		string dReactionName;
		int dAmplitudeCounter;
		ofstream dConfigFileStream;
		HCA_TensorContractor dHelicityCouplingAmplitudeBuilder;

		//Decay Trees
		NodeBase* gProductionNode; //top node

		//other
		set<Particle*> dFinalStateParticles;
		set<FitParamInfo*> dFitParams;
		map<Particle*, int> dParticleIndexMap; //int is prodcut index //same order as data
		map<int, Particle*> dReverseParticleIndexMap; //same as dParticleIndexMap, but key/value reversed
		map<vector<int>, pair<string, int> > gAmplitudeNames; //vector<int> are indices of independent amplitude terms, int is amp-counter name to constrain-to
		deque<string> dUniqueFullAmplitudeNames; //for initialization
		string dPIDString;

		void Generate_HelicityCouplingAmplitudes(deque<NodeBase*>& locTransitionTrees, deque<HCA_ContractionResult>& locAllFinalAmplitudes);
		inline ResonanceNode* Get_StartingNodes(NodeBase* locInputTransitionTree, deque<ResonanceNode*>& locOtherUndecayedNodesToProcess);

		void Generate_ConfigFile(deque<NodeBase*>& locTransitionTrees);
		void Build_HelicityPermutations(deque<deque<int> >& locHelicityTimes2Permutations);
		bool Handle_HelicityDecursion(int& locParticleIndex, deque<int>& locHelicityTimes2s, deque<int>& locResumeAtValues, deque<int>& locJTimes2s);

		bool Loop_PrintingAmplitudes(ResonanceNode* locCurrentNode, deque<ResonanceNode*> locOtherUndecayedNodesToProcess, string locSumName, deque<int>& locHelicityTimes2s, vector<string> locAmplitudeArgumentsVector, vector<string> locBWArgumentsVector, vector<int> locContractionIndices);
		void Calculate_PossibleHelicityRange(ResonanceNode* locCurrentNode, deque<int>& locHelicityTimes2s, int& locMinH12x2, int& locMaxH12x2);
		bool Check_IsNodeAmplitudeZero(ResonanceNode* locCurrentNode, int locH12x2, int locContractionIndex);

		string Get_ProductIndices(ResonanceNode* locIsParentNode);
		void Get_ProductIndices(ResonanceNode* locCurrentNode, vector<int>& locProductIndices);
};

inline void HCA_ConfigGenerator::Set_Particles(Particle* locPsiBeam, Particle* locPsiTarget, deque<Particle*>& locFinalStateParticles)
{
	dParticleIndexMap[locPsiBeam] = 0;
	dReverseParticleIndexMap[0] = locPsiBeam;
	dParticleIndexMap[locPsiTarget] = -1;
	dReverseParticleIndexMap[-1] = locPsiTarget;

	for(size_t loc_i = 0; loc_i < locFinalStateParticles.size(); ++loc_i)
	{
		Particle* locParticle = locFinalStateParticles[loc_i];
		dFinalStateParticles.insert(locParticle);
		dParticleIndexMap[locParticle] = loc_i + 1;
		dReverseParticleIndexMap[loc_i + 1] = locParticle;
	}
}

inline void HCA_ConfigGenerator::Set_FitParams(set<FitParamInfo*>& locFitParams)
{
	dFitParams = locFitParams;
}

inline string HCA_ConfigGenerator::Get_ProductIndices(ResonanceNode* locIsParentNode)
{
	vector<int> locProductIndices;
	Get_ProductIndices(locIsParentNode, locProductIndices);
	ostringstream locProductIndexStream;
	for(size_t loc_i = 0; loc_i < locProductIndices.size(); ++loc_i)
	{
		if(loc_i != 0)
			locProductIndexStream << "_";
		locProductIndexStream << locProductIndices[loc_i];
	}
	return locProductIndexStream.str();
}

inline void HCA_ConfigGenerator::Get_ProductIndices(ResonanceNode* locCurrentNode, vector<int>& locProductIndices)
{
	//traverse decays
	map<Particle*, int>::iterator locIterator = dParticleIndexMap.find(locCurrentNode->dNodeParticle);
	if(locIterator != dParticleIndexMap.end())
		locProductIndices.push_back(locIterator->second);
	else
	{
		Get_ProductIndices(locCurrentNode->dLeftFinalNode, locProductIndices);
		Get_ProductIndices(locCurrentNode->dRightFinalNode, locProductIndices);
	}
}

inline ResonanceNode* HCA_ConfigGenerator::Get_StartingNodes(NodeBase* locInputTransitionTree, deque<ResonanceNode*>& locOtherUndecayedNodesToProcess)
{
	//Initialize (Get starting node(s))
	ResonanceNode* locCurrentNode = dynamic_cast<ResonanceNode*>(locInputTransitionTree);
	if(locCurrentNode != NULL)
	{
		locOtherUndecayedNodesToProcess.push_back((ResonanceNode*)locCurrentNode->dPreviousDecayNode);
		return locCurrentNode;
	}

	//is exchange node
	ExchangeNode* locExchangeNode = dynamic_cast<ExchangeNode*>(locInputTransitionTree);
	Particle* locPsi1 = locExchangeNode->dLeftFinalNode->dNodeParticle;
	if(locPsi1->dIsResonanceFlag)
		locOtherUndecayedNodesToProcess.push_back(locExchangeNode->dLeftFinalNode);

	Particle* locPsi2 = locExchangeNode->dRightFinalNode->dNodeParticle;
	if(locPsi2->dIsResonanceFlag)
		locOtherUndecayedNodesToProcess.push_back(locExchangeNode->dRightFinalNode);

	if(locOtherUndecayedNodesToProcess.empty())
		return NULL;

	locCurrentNode = locOtherUndecayedNodesToProcess.back();
	locOtherUndecayedNodesToProcess.pop_back();
	return locCurrentNode;
}

#endif //_AmpToolsConfigGenerator_