#include "TFile.h"
#include "TTree.h"
#include "TTreeReader.h"
#include "TTreeReaderArray.h"
#include "TH2F.h"
#include "TStyle.h"
#include "TCanvas.h"
#include "TCutG.h"
#include "TChain.h"
#include <cstdlib>
#include <iostream>
#include <string>

#define NDET 16

TH1F * tDiff[NDET];
TH1F * tDiffg[NDET];
TH2F * PSD[NDET];
TH2F * PSDg[NDET];

TH2F * eL_tDiff[NDET];

ushort eL[NDET];
ushort eS[NDET];

const float timeResol = 2.0; //ns
const int timeRange[2] = {200, 800};
const int timeBin = (timeRange[1] - timeRange[0])/timeResol;

void analyzer(){

  for( int i = 0; i < NDET; i++ ){

    tDiff[i] = new TH1F(Form("Tdiff-%d", i), Form(" Detector %d", i), timeBin, timeRange[0], timeRange[1]);

    PSD[i] = new TH2F(Form("PSD-%d", i), Form(" Detector %d", i), 4096, 0, 30000,512,-1,1);

    eL_tDiff[i] = new TH2F(Form("E_tDiff-%d", i), Form(" Detector %d; TOF [ns]; eL", i), timeBin, timeRange[0], timeRange[1] ,1000, 0,15000);

  }

  TFile * file = new TFile("run268_3000.root");
  TTree * tree = (TTree *) file->Get("tree"); 

  TTreeReader reader(tree);

  TTreeReaderValue<ULong64_t>  evID = {reader, "evID"};
  TTreeReaderValue<UInt_t>    multi = {reader, "multi"};
  TTreeReaderArray<UShort_t>     sn = {reader, "sn"}; // serial no.
  TTreeReaderArray<UShort_t>     ch = {reader, "ch"}; // channel
  TTreeReaderArray<UShort_t>      e = {reader, "e"};  //long
  TTreeReaderArray<UShort_t>     e2 = {reader, "e2"}; //short
  TTreeReaderArray<ULong64_t>   e_t = {reader, "e_t"}; // in ns
  TTreeReaderArray<UShort_t>    e_f = {reader, "e_f"}; // in ps

  //^###########################################################
  //^ * Process
  //^###########################################################

  while(reader.Next()){

    unsigned long long tRF = 0;
    unsigned long long tN[NDET];
    for( int i =0; i < NDET; i++){
      tN[i] = 0;
      eL[i] = 0;
      eS[i] = 0;
    }

    if( *multi > 2) continue;

    for( int i = 0; i < *multi; i++){

      int bd = sn[i] - 6; // 6 = RF, 7 = Neutron detector
      int ID = ch[i];

      if( sn[i] == 6){
        tRF = e_t[i] * 1000 + e_f[i]; //in ps
      }else{
        continue;
      }

      eL[ID] = e[i];
      eS[ID] = e2[i];

      tN[ID] = e_t[i] * 1000 + e_f[i]; //in ps

    } // end multiplicity for loop

    for( int i = 0; i < NDET; i++) {
      if( tRF == 0 || tN[i] == 0 ) continue;

      double tdf;

      if( tRF > tN[i]){
        tdf = (tRF - tN[i])*1./1000.;
      }else{
        tdf = (tN[i] - tRF)*(-1.)/1000.;
      }

      //std::cout << "Finito" << std::endl;

      PSD[i]->Fill(eL[i],(eL[i]-eS[i])*1.0/eL[i]);
      tDiff[i]->Fill( tdf );

    } // end for loop over NDET

  }  // end while(reader.Next())



  gStyle->SetOptStat(11111111);

  TCanvas *canvas = new TCanvas("canvas", "anlayszer", 1000, 1000);
  canvas->Divide(4,4);

  for( int i = 0 ; i < NDET; i++ ){

    canvas->cd(i+1);
    tDiff[i]->Draw();

  }
  
} // end script()