#include "/home/astamato/Documents/lib_Thanos/header_root.h"
#include "/home/astamato/Documents/lib_Thanos/header_c.h"
#include <algorithm>


void add_coincidences(TString filename_in, TString filename_out){

//_____________________________________________________________________________________________________________________
//(0) Usefull variables

	// Proton, bunch counting
	bool new_T0 = false;

	
	// Misc
	double bin_content = 0.;//Will be used to fill the coincidence histogram, calculate time differences
	int    index_2     = 0;
	int    index_3     = 0;
	int    T0s         = 0;

	
//_____________________________________________________________________________________________________________________	
//(1) Open the root file and loop over the TTree

	TFile *f_root = TFile::Open(filename_in, "UPDATE");
	TTree *t      = (TTree*)f_root->Get("DICER");

	double         tof, tof0, dt;
	int            detn, RunNumber, event_ID;
	
	int   detn_prev, RunNumber_prev;
	float tof_prev, tof0_prev;
	
	vector<double> tof_2  , tof_3, DT;
	vector<int>    event_2, event_3, EVENT_ID, row_2, row_3, ROW;  


	t->SetBranchAddress("RunNumber_THANOS", &RunNumber);
	t->SetBranchAddress("event_ID_THANOS" , &event_ID);
	t->SetBranchAddress("tof0_THANOS"     , &tof0);
	t->SetBranchAddress("detn_THANOS"     , &detn);
	t->SetBranchAddress("tof_THANOS"      , &tof);
	
	//TBranch *bdt = t->Branch("dt", &dt);

	int     entries = t->GetEntries();
	
	std::cout << "Total entries found: " << entries << std::endl;
	
	float	j = entries/100.;
	int 	n = 1;
		
		
//_____________________________________________________________________________________________________________________		
//(3) Define the histograms
	
	// For the coincidneces
	TH1D *h_diff_2_3 = new TH1D("h_diff_2_3", "Difference 2 vs 3; Time (ns); Counts", 2000, -1000, 1000);
	TH1D *h_test     = new TH1D("h_test", "Test histogram", 2000, -1000, 1000);   

	n =1;
	printf("         10        20        30        40        50        60        70        80        90        100 %%\n");
	printf(" ---------|---------|---------|---------|---------|---------|---------|---------|---------|---------|    \n");
	
	
//_____________________________________________________________________________________________________________________
// () Loop over the entries
	
	for (int i = 0; i < entries; ++i){
		if(i == 0) {cout << " ";}
		if(i >= n*j) {cout << '|' << flush; n++;}// Progress indicator

		t->GetEntry(i);

		// All detectors apart from 2, 3 and 200
		if ( detn!=2 && detn!=3 && detn!=200 ){
			dt = -detn*1.e10;
			EVENT_ID.push_back(event_ID);
			DT.push_back(dt);
			ROW.push_back(i);
			continue;
		}
	
		// Detector 200 - check if there's a new T0
		if (detn==200){// new T0
		
			// Since we are in a new T0 let's calculate the differences in the previous one.
			if (T0s > 0){//make sure to ignore signals before the arrival of the 1st T0 
				if( event_2.size()>0 && event_3.size()>0 ){//make sure there are events in both detectors
					//Find the differences 2 - 3
					for (int ii=0; ii<event_2.size(); ++ii){
						for (int jj=0; jj<event_3.size(); ++jj){
							bin_content = (double) tof_2[ii] - tof_3[jj];//current difference
							h_diff_2_3->Fill(bin_content);
							
							if (jj==0){
								dt    = bin_content;// assign the minimum to the 1st difference
								index_3 = jj;
								index_2 = ii;
							}
							
							//Find differences 
							if ( abs(bin_content)<abs(dt) ){// |current| < |minimum|
								dt      = bin_content;
								index_3 = jj;
								index_2 = ii;
							}
						} 
						EVENT_ID.push_back( event_2[index_2] );
						DT.push_back( dt );	
						ROW.push_back( row_2[index_2] );	
					}
					// Now find the differences 3 - 2; Each signal must have a difference
					for (int jj=0; jj<event_3.size(); ++jj){
						for (int ii=0; ii<event_2.size(); ++ii){
							bin_content = (double) tof_3[jj] - tof_2[ii];//current difference
							
							if (ii==0){
								dt    = bin_content;// assign the minimum to the 1st difference
								index_2 = ii;
								index_3 = jj;
							}
							
							//Find differences 
							if ( abs(bin_content)<abs(dt) ){// |current| < |minimum|
								dt      = bin_content;
								index_3 = jj;
								index_2 = ii;
							}
						}
						EVENT_ID.push_back( event_3[index_3] );
						DT.push_back( dt );
						ROW.push_back( row_3[index_3] );
					}
				}
				
				if (event_2.size()==0 && event_3.size()>0){// Events only in 3
					dt = 3.e11;
					for (int jj=0; jj<event_3.size(); ++jj){
						EVENT_ID.push_back( event_3[jj] );
						DT.push_back( dt );
						ROW.push_back( row_3[jj] );
					}
				}
				
				if (event_2.size()>0 && event_3.size()==0){// Events only in 2
					dt = 2.e11;
					for (int ii=0; ii<event_2.size(); ++ii){
						EVENT_ID.push_back( event_2[ii] );
						DT.push_back( dt );
						ROW.push_back( row_2[ii] );
					}
				}
				
				
				event_2.clear();
				event_3.clear();
				tof_2.clear();
				tof_3.clear();
				row_2.clear();
				row_3.clear();
			
			}// T0s>0
		
			T0s++;
			dt = -detn*1.e10;
			EVENT_ID.push_back(event_ID);
			DT.push_back(dt);
			ROW.push_back(i);
			continue;
		}//detn==200
		
		if (detn == 2){
			if (T0s<=0){//Pulses before the 1st T0
				dt = detn*1.e10;
				EVENT_ID.push_back(event_ID);
				DT.push_back(dt);
				ROW.push_back(i);	
			}
			else{
				tof_2.push_back(tof);
				event_2.push_back(event_ID);
				row_2.push_back(i);
			}
			continue;
		}//detn==2
				
		if (detn==3){
			if (T0s<=0){//Pulses before the 1st T0
				dt = detn*1.e10;
				EVENT_ID.push_back(event_ID);
				DT.push_back(dt);
				ROW.push_back(i);
			}
			else{
				tof_3.push_back(tof);
				event_3.push_back(event_ID);
				row_3.push_back(i);
			}
			continue;
		}//detn==3
		
		
	}//_loop over entries
	
	//This is for the last T0
	int row_size = ROW.size();
	if (row_size < entries){
	
		if (event_2.size()==0 && event_3.size()>0){//Events only in 3
			dt = 3.e11;
			for (int i=0; i<event_3.size(); ++i){
				EVENT_ID.push_back( event_3[i] );
				DT.push_back( dt );
				ROW.push_back( row_3[i] );
			}
			event_3.clear();
			tof_3.clear();
			row_3.clear();
		}
		
		if (event_2.size()>0 && event_3.size()==0){// Events only in 2
			dt = 2.e11;
			for (int i=0; i<event_2.size(); ++i){
				EVENT_ID.push_back( event_2[i] );
				DT.push_back( dt );
				ROW.push_back( row_2[i] );
			}
			event_2.clear();
			tof_2.clear();
			row_2.clear();
		}
		
		if (event_2.size()>0 && event_3.size()>0){//Events in both detectors
			//Find the differences 2 - 3
					for (int ii=0; ii<event_2.size(); ++ii){
						for (int jj=0; jj<event_3.size(); ++jj){
							bin_content = (double) tof_2[ii] - tof_3[jj];//current difference
							h_diff_2_3->Fill(bin_content);
							
							if (jj==0){
								dt    = bin_content;// assign the minimum to the 1st difference
								index_3 = jj;
								index_2 = ii;
							}
							
							//Find differences 
							if ( abs(bin_content)<abs(dt) ){// |current| < |minimum|
								dt      = bin_content;
								index_3 = jj;
								index_2 = ii;
							}
						} 
						EVENT_ID.push_back( event_2[index_2] );
						DT.push_back( dt );	
						ROW.push_back( row_2[index_2] );	
					}
					// Now find the differences 3 - 2; Each signal must have a difference
					for (int jj=0; jj<event_3.size(); ++jj){
						for (int ii=0; ii<event_2.size(); ++ii){
							bin_content = (double) tof_3[jj] - tof_2[ii];//current difference
							
							if (ii==0){
								dt    = bin_content;// assign the minimum to the 1st difference
								index_2 = ii;
								index_3 = jj;
								//place = 
							}
							
							//Find differences 
							if ( abs(bin_content)<abs(dt) ){// |current| < |minimum|
								dt      = bin_content;
								index_3 = jj;
								index_2 = ii;
							}
						}
						EVENT_ID.push_back( event_3[index_3] );
						DT.push_back( dt );
						ROW.push_back( row_3[index_3] );
					}
					event_2.clear();
					event_3.clear();
					tof_2.clear();
					tof_3.clear();
					row_2.clear();
					row_3.clear();
		}
		
		std::cout << std::endl;

	}
	
	std::cout << std::endl;
	std::cout << "Sorting... " << std::endl;
	std::sort( ROW.begin(), ROW.end() );
	std::cout << "Finished sorting!" << std::endl;
	std::cout << std::endl;
	
	
	for (int i = 0; i < entries; ++i){
		if(i == 0) {cout << " ";}
		if(i >= n*j) {cout << '|' << flush; n++;}// Progress indicator

		t->GetEntry(i);
		// Now I have somehow to sort the coincidences with respect to THE ROW.
		if (i != ROW[i]){
			std::cout << "row = " << i << std::endl;
			break;
		}
	}
	

	h_diff_2_3->Write();
	h_test->Write();
	f_root->Close();
	
	std::cout << "Entries = " << entries << std::endl;
	std::cout << "Size    = " << EVENT_ID.size() << std::endl;
	std::cout << "Rows    = " << ROW.size() << std::endl;
	EVENT_ID.clear();
	DT.clear();
	ROW.clear();
	std::cout << "Size    = " << EVENT_ID.size() << std::endl;

}