#include "TFile.h"
#include "TH1.h"
#include "TH2.h"
#include "TChain.h"
#include "TTree.h"
#include "TLeaf.h"
#include <ROOT/RDataFrame.hxx>

#include <fstream>
#include <stdio.h>
#include <iostream>
#include <vector>
#include <bitset>

#define UpdateEvents    20000

#ifndef CALIBRATION
// default setting:---> 0.5 keV/ch
#define CALIBRATION 0.5
#endif

#if !defined(No_Clover)
// default setting: No_Clover = 32
#define No_Clover 32
#endif /* !defined(No_Clover) */

using namespace std;

bool isPositive(double x){return x > 0;}

TFile *fop = NULL;

TH1D *his[150];
TH1D *his_cal[150];
TH1D *his_AB[32];
TH1D *his_TDC[32];
TH1D *his_TAC;

float coffset[150] = {0.0}, tsqrt[150] = {0.0}, slope[150] = {0.5}, quad[150] = {0.0}, nquad[150] = {0.0};
signed long Eclab[No_Clover], CLT[No_Clover], energy[150], energy_cal[150];

int Did[No_Clover], nOFlowTDC=0;
unsigned int TLocation[No_Clover];
signed int TOffset[No_Clover] = {0};
int cryst_index[No_Clover][4] = {9999};

vector<unsigned long long> vCLT;

ifstream infile("Names.dat");		// File containing histogram names, one per line
ifstream DetFile2("detector.dat");	// Energy & Time signal mapping for each detector 
ifstream calfile2("calib.dat");		// File with calibration parameters

//std::vector<ROOT::RDF::RResultPtr<TH1D>> GetValue(const char *filename = "eu152_11july_afterGlitch_002.root", const char *treename = "RoseNIAS")
std::vector<ROOT::RDF::RResultPtr<TH1D>> GetValue(const char *filename = "eu152_11july_afterGlitch_*0*.root", const char *opfname = "test.root")
//void GetValue(const char *filename = "eu152_11july_afterGlitch_*0*.root", const char *treename = "RoseNIAS")
//void GetValue(const char *filename = "eu152_11july_afterGlitch_*0*.root", const char *opfname = "test.root")
{
ROOT::EnableImplicitMT();

	fop = new TFile(opfname, "RECREATE", "", 101);

	const char *treename = "RoseNIAS";

	Int_t NoPara = 105;		//No of parameters in experiment. Also, number of branches in Tree "RoseNIAS"
	char nm[30];
	const char *nmHis[20];
	string nmRawHis[105];
	string nmCalHis[105];

/*
	for (int i=0; i<NoPara; ++i){
		infile>>nm;
		his[i]=new TH1D(nm,nm,8192,0.,8192.);		//Raw histograms
		nmRawHis[i] = nm;
	}
	infile.seekg(0, infile.beg);

	for (int i=0; i<NoPara; ++i){
		infile>>nm;
		strcat(nm,"_cal");
		his_cal[i]=new TH1D(nm,nm,8192,0.,8192.);	//Calibrated histograms 
		nmCalHis[i] = nm; 
	}

	for (int i = 0; i < No_Clover; i++){
        snprintf(nm, sizeof(nm), "TDCA_%02i", i+1);		//Aligned TDCs
        his_TDC[i]=new TH1D(nm,nm,8192,0.,8192.);
	}

	for (int i = 0; i < No_Clover; i++){
        snprintf(nm, sizeof(nm), "AB%02i", i+1);		//Add-back clover
        his_AB[i]=new TH1D(nm,nm,8192,0.,8192.);
	}
	
	his_TAC=new TH1D("TACspectra","TACspectra",8192,0.,8192.);	//TAC Spectrum
*/

	for(Int_t l = 0; l < NoPara; l++){
	//Reading Calibration Parameters (obtained from 2nd order polynomial fit)
        //calfile2 >> coffset[l] >>  slope[l] >> quad[l];
	//Reading Calibration Parameters (obtained from 3rd order polynomial fit)
        //calfile2 >> coffset[l] >>  slope[l] >> quad[l] >> nquad[l];
	//Reading Calibration Parameters (obtained from 2nd order polynomial + sqrt term fit)	
        calfile2 >> coffset[l] >> tsqrt[l] >> slope[l] >> quad[l];	
	}

	//Reading detector mapping
	int DcrystIndex, DTLo, DTOff;
	for(int l = 0; l < No_Clover; l++){
		DetFile2 >> Did[l];
		for(int k = 0; k < 4; k++){
                DetFile2 >> DcrystIndex;
                cryst_index[l][k] = DcrystIndex -1;
		}
		DetFile2 >> DTLo >> DTOff;
		TLocation[l] = DTLo - 1;  TOffset[l] = DTOff;
		//cout << "Clover # [" << Did[l] << "]: \t ADC Mapping \t" << cryst_index[l][0] + 1 << " \t " << 
		//cryst_index[l][1] + 1 << "\t" << cryst_index[l][2] + 1 << "\t" << cryst_index[l][3] + 1 << 
		//" \t TDC Location:   " << DTLo << " \t TDC Offset:   " << DTOff << endl;
	}
	
	ROOT::RDataFrame df(treename, filename);

	long unsigned nEvents = df.Count().GetValue();
	cout << "Number of events to be processed: \t" << nEvents << endl;

	std::vector<ROOT::RDF::RResultPtr<TH1D>> histos;
	for (const std::string &col : df.GetColumnNames()) {
		ROOT::RDF::TH1DModel hmod1D(col.c_str(), col.c_str(), 8192, 0, 8192);
		auto h = df.Filter(isPositive, {col}).Histo1D<double>(hmod1D, col);
		histos.push_back(h);
		auto vectorPtr = df.Take<double>(col);
	}
	auto eventCount = df.Count();
        eventCount.OnPartialResult(10000,
                                [](ULong64_t e) {std::cout << e << " : events processed " <<  std::flush << endl; });


	histos[37]->Draw();
	histos[37]->GetMean();


fop->cd();
fop->Write();
//delete fop; // automatically deletes the tree and all histograms

return histos;

}