/*
  https://root.cern.ch/doc/master/classTAttLine.html
 */
// compile with:
// g++ -o plotevent-mcmx.exe plotevent-mcmx.cpp `root-config --cflags --glibs`
//
#include "TH1F.h"
#include "TH2F.h"
#include "TH3F.h"
#include "TProfile.h"
#include "TMath.h"
#include "TF1.h"
#include "TLegend.h"
#include "TCanvas.h"
 //
#include "TROOT.h"
#include "TStyle.h"
#include "TTree.h"
#include "TFile.h"
#include "TGraphErrors.h"
#include "TStyle.h"
#include <string>
#include <cmath>
#include <iostream>
#include <fstream>
#include <vector>
#include <algorithm>
//
#include <sys/types.h>
#include <sys/stat.h>
#include <unistd.h>

using namespace std;


int plotevent(int thr = 0,int event = 1, const char *dir = "",  const char *filename = "", const char *title = "") {
  char name[100];
  char buffer[512];
  Char_t iChan; // Channel number
  Char_t iChip; // Channel number
  UShort_t iMcm; // MCM/FEC number
  UShort_t iClLength; //number of samples
  UShort_t iClTime;
  UInt_t iEventnb;
  Int_t iSamples[1024]; // Samples
  int nChan = 0;
  Int_t mappingx[256],mappingy[256];
  TH1F* hPed;
  TH2F* hPed2;
  TH1F* hAdc;
  TH2F* hHit;
  TH2F* hTimeProj;
  TH3F* hTime;
  TH1F* hTimeHits;
  TH1F* hEventnbsweep;
  TH1F* hTriggersweep;
  TH1F* hTimesweeps[384];
  int lastSamples = 20;
  
  // read mapping file
  FILE *mapfile;
  sprintf(name,"%s","mapping-2.csv");
  mapfile = fopen(name,"r");
  if (mapfile == NULL) {
    printf("Mapfile %s not found\n",name);
      return 1;
  }
  int line = 0;
  int row,column;
  int chip,channel;
  int mapindex;
  int dcl;
  char *ptr;
  memset(mappingx,-1,sizeof(mappingx));
  memset(mappingy,-1,sizeof(mappingy));
  while (fgets(buffer,sizeof(buffer),mapfile) != NULL) {
    if (line != 0) {
      //      puts(buffer);
      ptr = strtok(buffer,";\n");
      column=0;
      while(ptr!=NULL) {
	if (column == 0) sscanf(ptr,"%d",&dcl);
	else if (column == 1) sscanf(ptr,"%d",&row);
	else {
	  if (strcmp(ptr,"GND") != 0) {
	    sscanf(ptr,"%d,%d",&chip,&channel);
	    mapindex = dcl*128 + chip*16 + channel;
	    if ((mappingx[mapindex] >= 0) || (mappingy[mapindex] >= 0)) {
	      printf("MAP %d(%d,%d) row %d col %d: %d %d\n",mapindex,chip,channel,row,column-2,mappingx[mapindex],mappingy[mapindex]);
	    }
	    else {
	      mappingx[mapindex] = column - 2;
	      mappingy[mapindex] = row;
	    }
	    //	    printf("MAP %d %d %d\n",mapindex,row,column-2);
	  }
	  else printf("GND %d %d\n",row,column-2);
	}
	ptr = strtok(NULL,";\n");
	column++;
      }
    }
    line++;
  }
  fclose(mapfile);
  mapfile = fopen("mapping.txt","w");
  for (int i = 0; i < 256; i++) {
    fprintf(mapfile,"INDEX: %d: x:%d y:%d\n",i,mappingx[i],mappingy[i]);
  }
  fclose(mapfile);
    /* 
  for (int j = 0; j < 256; j++) {
    printf("%d: X(col) Y(row) %d %d\n",j,mappingx[j],mappingy[j]);
  }
  */
  hPed = ((TH1F *)(gROOT->FindObject("hPed")));
  if (hPed) delete hPed;
  hPed = new TH1F("hPed", "Pedsub",40,-20.0,20.0);
  hPed->SetMinimum(0);
  hPed->GetXaxis()->SetTitle("Pedsub");
  hPed->GetYaxis()->SetTitle("Entries");

  hPed2 = new TH2F("hPed2", "Pedsub",384,0,384,100,50,150);
  hPed2->SetMinimum(0);
  hPed2->GetXaxis()->SetTitle("Pedestal");
  hPed2->GetYaxis()->SetTitle("Channel");

  hAdc = ((TH1F *)(gROOT->FindObject("hAdc")));
  if (hAdc) delete hAdc;
  hAdc = new TH1F("hAdc", "Adc",1024,0.0,1024.0);
  hAdc->SetMinimum(0);
  hAdc->GetXaxis()->SetTitle("ADC");
  hAdc->GetYaxis()->SetTitle("Entries");
 
  
  hHit = ((TH2F *)(gROOT->FindObject("hHit")));
  if (hHit) delete hHit;
  sprintf(name,"XY EVENT %d",event);
  hHit = new TH2F("hHit",name,26,-0.5,25.5,10,-0.5,9.5);
  hHit->SetMinimum(0);
  hHit->SetMaximum(500);
  hHit->GetXaxis()->SetTitle("Column");
  hHit->GetYaxis()->SetTitle("Row");

  hTimeHits = new TH1F("hTimeHits", "Time hits",1024,-0.5,1023.5);
  hTimeHits->SetMinimum(0);
  hTimeHits->GetXaxis()->SetTitle("Time");
  hTimeHits->GetYaxis()->SetTitle("Entries");

  for (int i = 0; i < 384; i++) {
    sprintf(name,"hTimesweeps_%d",i);
    hTimesweeps[i] = new TH1F(name, "Time sweep",1024,-0.5,1023.5);
    sprintf(name,"MCM %d CHIP %d chan %d",i/128,(i/16)%8,i%16);
    hTimesweeps[i]->SetTitle(name);
  }

  // Event number timesweep
  int EventnbChan = 128*2 + 16*3;
  hEventnbsweep = new TH1F("hEventnbsweep", "Event  number sweep",1024,-0.5,1023.5);
  hEventnbsweep->SetTitle("Event number sweep");
  // Trigger timesweep
  int TriggerChan = 128*2 + 16*3 + 15;
  hTriggersweep = new TH1F("hTriggersweep", "Trigger time sweep",1024,-0.5,1023.5);
  hTriggersweep->SetTitle("Trigger time sweep");
  
  
  hTimeProj = ((TH2F *)(gROOT->FindObject("hTimeProj")));
  if (hTimeProj) delete hTimeProj;
  sprintf(name,"TimeProjection EVENT %d",event);
  hTimeProj = new TH2F("hTimeProj",name,1024,-0.5,1023.5,10,-0.5,9.5);
  hTimeProj->GetXaxis()->SetTitle("Time");
  hTimeProj->GetYaxis()->SetTitle("Row");
  

  hTime = ((TH3F *)(gROOT->FindObject("hTime")));
  if (hTime) delete hTime;
  sprintf(name,"Time View EVENT %d",event);
  //  hTime = new TH3F("hTime", "Time",26,-0.5,25.5,10,-0.5,9.5,200,0.0,1024.0);
  hTime = new TH3F("hTime",name,1024,0.0,1024.0,26,-25.5,0.5,10,-0.5,9.5);
  hTime->SetMinimum(0);
  hTime->SetMaximum(1.1);
  hTime->GetXaxis()->SetTitle("Time sample");
  hTime->GetXaxis()->SetTitleOffset(1.5); // offset = 1.5 will add 50 per cent more to the standard offset.
  hTime->GetYaxis()->SetTitle("Column");
  hTime->GetZaxis()->SetTitle("Row");
  hTime->GetYaxis()->SetBinLabel(26,"0");
  hTime->GetYaxis()->SetBinLabel(21,"5");
  hTime->GetYaxis()->SetBinLabel(16,"10");
  hTime->GetYaxis()->SetBinLabel(11,"15");
  hTime->GetYaxis()->SetBinLabel(6,"20");
  hTime->GetYaxis()->SetBinLabel(1,"25");
  hTime->GetZaxis()->SetBinLabel(10,"0");
  hTime->GetZaxis()->SetBinLabel(9,"1");
  hTime->GetZaxis()->SetBinLabel(8,"2");
  hTime->GetZaxis()->SetBinLabel(7,"3");
  hTime->GetZaxis()->SetBinLabel(6,"4");
  hTime->GetZaxis()->SetBinLabel(5,"5");
  hTime->GetZaxis()->SetBinLabel(4,"6");
  hTime->GetZaxis()->SetBinLabel(3,"7");
  hTime->GetZaxis()->SetBinLabel(2,"8");
  hTime->GetZaxis()->SetBinLabel(1,"9");
  hHit->GetYaxis()->SetBinLabel(10,"0");
  hHit->GetYaxis()->SetBinLabel(9,"1");
  hHit->GetYaxis()->SetBinLabel(8,"2");
  hHit->GetYaxis()->SetBinLabel(7,"3");
  hHit->GetYaxis()->SetBinLabel(6,"4");
  hHit->GetYaxis()->SetBinLabel(5,"5");
  hHit->GetYaxis()->SetBinLabel(4,"6");
  hHit->GetYaxis()->SetBinLabel(3,"7");
  hHit->GetYaxis()->SetBinLabel(2,"8");
  hHit->GetYaxis()->SetBinLabel(1,"9");
  
  hTimeProj->GetYaxis()->SetBinLabel(10,"0");
  hTimeProj->GetYaxis()->SetBinLabel(9,"1");
  hTimeProj->GetYaxis()->SetBinLabel(8,"2");
  hTimeProj->GetYaxis()->SetBinLabel(7,"3");
  hTimeProj->GetYaxis()->SetBinLabel(6,"4");
  hTimeProj->GetYaxis()->SetBinLabel(5,"5");
  hTimeProj->GetYaxis()->SetBinLabel(4,"6");
  hTimeProj->GetYaxis()->SetBinLabel(3,"7");
  hTimeProj->GetYaxis()->SetBinLabel(2,"8");
  hTimeProj->GetYaxis()->SetBinLabel(1,"9");

  TCanvas *c1 = new TCanvas("c1","Single cluster",600,600);
 
  sprintf(name, "%s/%s.root",dir,filename);
  printf("Run: %s\n",name);

  TFile *inFile = new TFile(name);
  TTree *tree = (TTree*) inFile->Get("tree");
  if(!tree) {
    cerr << "Tree SALTRO events was not found: " << name << "\n";
    return 1;
  }
  // Assign variables to branches
  tree->SetBranchAddress("iEventnb",&iEventnb);
  tree->SetBranchAddress("iSamples",&iSamples);
  tree->SetBranchAddress("iChan",&iChan);
  tree->SetBranchAddress("iChip",&iChip);
  tree->SetBranchAddress("iMcm",&iMcm);
  tree->SetBranchAddress("iClLength",&iClLength);
  tree->SetBranchAddress("iClTime",&iClTime);

  Int_t nentries = tree->GetEntries(); 
  printf("ENTRIES %d\n",nentries);
  Int_t nchan;
  sprintf(name,"%s/%s-pedfile.csv",dir,filename);
  FILE *csvpedfile = fopen(name,"w");
  if (event == 0) {
    fprintf(csvpedfile,"iEventnb,nChan,ped\n");
  }
  for (Int_t j = 0; j < nentries; j++) {
    tree->GetEntry(j);
    if ((event == iEventnb) || (event == 0)) {
      //      printf("EVENTNB %d %d %d %d\n",iEventnb,iMcm,iChip,iChan);
      nChan = iMcm*128 + iChip * 16 + iChan;
      if (iClLength > 0) {
	float ped = 0;
	int pedsamples = 0;
	if (thr > 0) {
	  //	  for (int i = (iClLength-lastSamples); i < (iClLength-1); i++) {
	  for (int i = 0; i < lastSamples; i++) {
	    ped += iSamples[i];
	    pedsamples++;
	  }
	  ped = ped/float(pedsamples);
	  if (event == 0) {
	    fprintf(csvpedfile,"%d,%d,%.3f\n",iEventnb,nChan,ped);
	  }
	  hPed2->Fill(nChan,ped);
	}
	if (iMcm <= 1) {
	  row = mappingy[nChan];
	  column = mappingx[nChan];
	  int maxval = 0;
	  //printf("LENGTH1 %d TIME %d\n",iClLength,iClTime);
	  for (int i = 0; i < iClLength-1; i++) {
	    hAdc->Fill(iSamples[i]);
	    iSamples[i] -= (int)ped;
	    //	  if ( i >= (iClLength-lastSamples)) hPed->Fill(iSamples[i]);
	    if ( i < lastSamples) hPed->Fill(iSamples[i]);
	    if (iSamples[i] < 0) iSamples[i] = 0;
	    //	    if (iSamples[i] > thr) {
	      hHit->Fill(column,9-row);
	      int chan = iMcm*128 + iChip*16 + iChan;
	      if ((iClTime-i) < 0) {
		printf("EVENTNB %d LENGTH %d TIME %d INDEX %d INDEX-i %d VALUE %d CHAN %d M%d C%d CH%d\n",iEventnb,iClLength,iClTime,i,iClTime-i,iSamples[i],chan, iMcm,iChip,iChan);
	      }
	      else {
		hTimesweeps[chan]->AddBinContent(iClTime-i,iSamples[i]);
		//printf("LENGTH %d %d %d\n",i,iSamples[i],chan);
		hTimeHits->Fill(iClTime-i);
		hTimeProj->Fill(iClTime-i,9-row,float(iSamples[i])/1024.0);
		hTime->Fill(iClTime-i,-column,9-row,float(iSamples[i])/1024.0);
	      }
	      //	      if (event == 0) printf("EVT %d nChan %d sample %d value %d\n",iEventnb,nChan,i,iSamples[i]);
	      //}
	  }
	  //printf("LENGTH2 %d\n",iClLength);
	}
	else if (iMcm ==2) {
	  for (int i = 0; i < iClLength-1; i++) {
	    int chan = iMcm*128 + iChip*16 + iChan; 
	    //	    if (iSamples[i] < 0) iSamples[i] = 0;
	    //	    if (iSamples[i] > thr) {
	      hTimesweeps[chan]->AddBinContent(iClTime-i,iSamples[i]);
	      if (chan == EventnbChan) {
		hEventnbsweep->AddBinContent(iClTime-i,iSamples[i]);
	      }
	      else if (chan == TriggerChan) {
		hTriggersweep->AddBinContent(iClTime-i,iSamples[i]);
	      }
	      //}
	  }
	}
      }
    }
  }
  printf("DONE\n");
  gStyle->SetOptStat(0);
  sprintf(name,"%s-evt%d-%s.pdf(",filename,event,title);
  puts(name);
  hHit->Draw("COLZ");
  c1->SaveAs(name,"pdf");
  /*
Float_t theta = 45;
Float_t phi   = 0;
gPad->SetPhi(phi);
gPad->SetTheta(theta);
  */
  hTimeProj->Draw("BOX");
  //  hTimeProj->SetAxisRange(0, 550, "X");
    hTimeProj->SetAxisRange(0, 1024, "X");
  hTimeProj->Draw("BOX");
  sprintf(name,"%s-evt%d-%s.pdf",filename,event,title);
  c1->SaveAs(name,"pdf");
  hTime->Draw("BOX");
  //hTime->SetAxisRange(0, 550, "X");
  hTime->SetAxisRange(0, 1024, "X");
  sprintf(name,"%s-evt%d-%s.pdf",filename,event,title);
  c1->SaveAs(name,"pdf");

  hTimeHits->Draw();
  c1->SaveAs(name,"pdf");
  hPed2->Draw();
  c1->SaveAs(name,"pdf");

  if (event == 0) hEventnbsweep->Scale(1.0/iEventnb);
  hEventnbsweep->Draw();
  c1->SaveAs(name,"pdf");
  if (event == 0) hTriggersweep->Scale(1.0/iEventnb);
  hTriggersweep->Draw();
  c1->SaveAs(name,"pdf");

  gStyle->SetOptStat(111110);
  hAdc->SetAxisRange(0, 4000, "Y");
  hAdc->SetAxisRange(0, 1024, "X");
  hAdc->Draw();
  sprintf(name,"%s-evt%d-%s.pdf)",filename,event,title);
  c1->SaveAs(name,"pdf");
  
  gStyle->SetOptStat(111111);
  hPed->Draw();
  sprintf(name,"%s-ped-evt%d-%s.pdf",filename,event,title);  
  c1->SaveAs(name,"pdf");
  fclose(csvpedfile);

  sprintf(name,"%s-ped-evt%d-%s.pdf",filename,event,title);  
  c1->SaveAs(name,"pdf");
  c1->Clear();
  c1->Divide(4,4);
  gStyle->SetOptStat(0);
  for (int i = 0; i < 24; i++) {
    int j;
    for (j = 0; j < 16; j++) {
      c1->cd(j+1);
      int chan = i*16 + j;
      //      hTimesweeps[chan]->SetAxisRange(60,120, "Y");
      if (event == 0) hTimesweeps[chan]->Scale(1.0/iEventnb);
      hTimesweeps[chan]->Draw();
    }
    if (i==0)   sprintf(name,"timesweeps.pdf(");
    else if (i==23)  sprintf(name,"timesweeps.pdf)");
    else  sprintf(name,"timesweeps.pdf");
    c1->SaveAs(name,"pdf");
  }


  return 0;
}

int main (int argc, char *argv[])
{
  int iret = 0;
  printf("%s %s %s %s %s\n",argv[1],argv[2],argv[3],argv[4],argv[5]);
  if (argc != 6) {
    printf("Wrong number of arguments: argc = %d\nSyntax should be:\n./plotevtmcm.exe event dir filename title\n",argc);
    return -1;
  }
  int thr = atoi(argv[1]);
  int event = atoi(argv[2]);
  iret =plotevent(thr,event,argv[3],argv[4],argv[5]);
  return iret;
}