#include <stdio.h>
#include <stdlib.h>
#include <iostream>
#include "vector"

//  root includes
#include "TFile.h"
#include "TTree.h"
#include "TCanvas.h"
#include "TH1D.h"
#include "TH2D.h"
#include "TGraph.h"
#include <TChain.h>

// decoding
#include "decoding.h"
#include "analysis.h"





//TTree* readTTree(const Char_t *inputdir, const Char_t *filename);

void ana_data(char* inputdir, char* outputdir, int pede, int check_max_subrun){
   Double_t **pedestal_X = pedestal_calc(inputdir,outputdir,pede, check_max_subrun, "X");
     
}


void run(){
   char inputdir[]= "/local/home/fbelloni/WORK/IPM/DATA/STRIPS/5_12_2017";                  // Dir where the TTrees are stored
   char outputdir[]="/local/home/fbelloni/WORK/IPM/DATA/STRIPS/Routines/Process/Complete/";  // Change later
   
   // /home/irfulx157/mnt/fbelloni/workspace/XYMG/ROUTINES/";    // dir where are stored the TTrees
   //int run = 15168;                        // run number    13833  13848
   int pede = 1;                       // pedestal run of reference
   int check_max_subrun = 1;
   //int lower_segment = 0;                    // lower segment
   //int upper_segment = 0;                    // higher segment
   //int nsegm = 0;
   //int stream_number = 2;                    // stream number containing detector
   //char MyDetName[] = "XYMG";                // detector name, detectors are FIC1, SILI, MGAS, XYMG, PKUP, ...
   //char Run_Type[] = "NORM";
   //Option_t* signals_status = "NSAT";    //SATU or NSAT (saturated or non saturated)
   ana_data(inputdir, outputdir, pede, check_max_subrun);//, run, pede, check_max_segm, signals_status);

}







int main(int argc, char **argv) {
   Bool_t dobatch = kFALSE; // kTRUE for batch without windows, kFALSE with windows
   if (dobatch) {
      gROOT->SetBatch();          // no window display, for batch processing, use TApplication instead of TRint
      TApplication *theApp = new TApplication("App", &argc, argv);
      run();
      return 0;
   } else {
      TRint *theApp = new TRint("App", &argc, argv, NULL, 0);
      run();
      theApp->Run();
      return 0;
   }
}



/*
void ProcessAll(){
    
 // Create a Canvas to "draw" on it   
 TCanvas*c1 = new TCanvas("c1","Histogram of strip pedestal signal",800,800);// 1900,1200);
 c1->Divide(3,2);
 TCanvas*c2 = new TCanvas("c2","Histogram of strip signal",800,800);// 1900,1200);
 c2->Divide(3,2);
  
 TCanvas*c3 = new TCanvas("c3","Histogram of strip signal - pedestal signal", 800,800);
 c3->Divide(3,2);

  TH1::AddDirectory(kFALSE);
  TH2::AddDirectory(kFALSE);

 // The following vextor is  used for defining the x axis in graphs
  double *xaxis= new double[32] {0.4 , 1.32 , 2.24 , 3.16 , 4.08 , 5 , 5.92 , 6.84 , 7.76 , 8.68 , 9.6 , 10.52 , 11.44 , 12.36 , 13.28 , 14.2 , 
                             15.12 , 16.04 , 16.96 , 17.88 , 18.8 , 19.72 , 20.64 , 21.56 , 22.48 , 23.4 , 24.32 , 25.24 , 26.16 , 27.08 , 28 , 28.92 };

 //===========================================
 //
 //      Treat the pedestal file 
 //                             
 //===========================================
                             
                             
 TFile*f_p = new TFile("../../5_12_2017/Pedestal.root","read");       // Open the root file containing the pedestal data in the form of a TTrre
 TTree*t_p = (TTree*)f_p->Get("DataTree");                            // Create a pointer to the TTree
 
 vector<double> *q_p           =0;                                    // Create a vector to contain the 32 signals for each one of the 32 strips
 vector<int>    *sat_p         =0;                                    // Create a vector to contain the 32 indexes saying, for each of the 32 strips, if the signal is saturated
 double          time_ns_p     =0;                                    // Variable to read the time stamp
 int             card_nb_p     =0;                                    // Variable to read the card number
   
 TBranch *charge_p     = 0;                                           // Create pointers to the branches of the TTree
 TBranch *saturation_p = 0;
 TBranch *card_p       = 0;
 TBranch *time_p       = 0;
 
 t_p->SetBranchAddress("time",&time_ns_p);                           // Get the branches
 t_p->SetBranchAddress("card",&card_nb_p);
 t_p->SetBranchAddress("Q",&q_p,&charge_p);
 t_p->SetBranchAddress("S",&sat_p,&saturation_p);
 
 TH1D*hX_p_nd = new TH1D("hX_p_nd","Strips X - decoding not applied", 32,0,32);                  // Create hitsograms of hX_p_nd = X strips, pedestal signal, not decoded (= decoding not applied)
 TH1D*hY_p_nd = new TH1D("hY_p_nd","Strips Y - decoding not applied", 32,0,32);                  // Create hitsograms of hY_p_nd = Y strips, pedestal signal, not decoded (= decoding not applied)

 TH1D*hX_p = new TH1D("hX_p","Strips X - decoding applied", 32,0,32);                    // Create hitsograms of hX_p = X strips, pedestal signal 
 int EntriesCard1_p = 0;
 TH1D*hY_p = new TH1D("hY_p","Strips Y - decoding applied", 32,0,32);                    // Create hitsograms of hY_p = Y strips, pedestal signal 
 int EntriesCard2_p = 0;
 
 
 

 int nEntries_p=t_p->GetEntries(); 

 for (Int_t i = 0; i < nEntries_p; i++) {                             // This loop is used to read the TTree data and fill the histograms. 
      
      t_p->GetEntry(i);

      //cout << "time " << time_ns_p << endl;
      //cout << "card_nb " << card_nb_p << endl;


      Long64_t tentry_p = t_p->LoadTree(i);
      charge_p->GetEntry(tentry_p);
      saturation_p->GetEntry(tentry_p);
      //for (UInt_t j = 0; j < q->size(); ++j) {
      // printf("%lf   %d \n",  q->at(j) , sat->at(j));}
      if(card_nb_p==1)   {EntriesCard1_p++;
                          for (UInt_t j = 0; j < q_p->size(); ++j) {
                                              hX_p_nd->Fill(j, q_p->at(j));
                                              hX_p->Fill(strip_xdec[j]-1,q_p->at(j));    //decoding applied here
                                                                   }
                        }
 
      if(card_nb_p==2)   {EntriesCard2_p++;
                          for (UInt_t j = 0; j < q_p->size(); ++j) {
                                              hY_p_nd->Fill(j, q_p->at(j));
                                              hY_p->Fill(strip_ydec[j]-1,q_p->at(j));}   // decoding applied here
                                              
                       }

   }

 
  // Plots
  c1->cd(1);
  double MaxXCont_p = 0;
  for(int i = 0; i<32; i++) {if(hX_p_nd->GetBinContent(i)>MaxXCont_p) MaxXCont_p = hX_p_nd->GetBinContent(i);}
  hX_p_nd->GetXaxis()->SetTitle("Strip number");
  hX_p_nd->GetYaxis()->SetTitleOffset(1.5);
  hX_p_nd->GetYaxis()->SetTitle("Q (nC)");
  hX_p_nd->SetLineColor(2);
  hX_p_nd->SetFillColor(2);
  hX_p_nd->GetYaxis()->SetRangeUser(0,MaxXCont_p+MaxXCont_p*0.4);
  hX_p_nd->SetBarWidth(0.95);
  hX_p_nd->Draw("BAR");
  cout<<"Events in card 1 pedestal " << EntriesCard1_p << endl;
  c1->cd(2);
  hX_p->GetXaxis()->SetTitle("Strip number");
  hX_p->GetYaxis()->SetTitleOffset(1.5);
  hX_p->GetYaxis()->SetTitle("Q (nC)");
  hX_p->SetLineColor(2);
  hX_p->SetFillColor(2);
  hX_p->GetYaxis()->SetRangeUser(0,MaxXCont_p+MaxXCont_p*0.4);
  hX_p->SetBarWidth(0.95);
  hX_p->Draw("BAR");
  c1->cd(3);
  double *x_graph_p_y = new double[32];
  for(int i=0;i<32;i++) x_graph_p_y[i] =hX_p->GetBinContent(i+1);  
  TGraph* gX_p_mm = new TGraph(32, xaxis, x_graph_p_y);
  gX_p_mm->SetTitle("Strips X - decoding applied");
  gX_p_mm->GetXaxis()->SetTitle("Detector side (mm)");
  gX_p_mm->GetYaxis()->SetTitle("Q (nC)");
  gX_p_mm->GetYaxis()->SetTitleOffset(1.5);
  gX_p_mm->GetYaxis()->SetRangeUser(0,MaxXCont_p+MaxXCont_p*0.4);
  gX_p_mm->SetLineColor(2);
  gX_p_mm->SetMarkerStyle(20);
  gX_p_mm->SetMarkerSize(1);
  gX_p_mm->Draw("APL");
  
  
  
  c1->cd(4);
  double MaxYCont_p = 0;
  for(int i = 0; i<32; i++) {if(hY_p_nd->GetBinContent(i)>MaxYCont_p) MaxYCont_p = hY_p_nd->GetBinContent(i);}
  hY_p_nd->GetXaxis()->SetTitle("Strip number");
  hY_p_nd->GetYaxis()->SetTitleOffset(1.5);
  hY_p_nd->GetYaxis()->SetTitle("Q (nC)");
  hY_p_nd->SetLineColor(2);
  hY_p_nd->SetFillColor(2);
  hY_p_nd->GetYaxis()->SetRangeUser(0,MaxYCont_p+MaxYCont_p*0.4);
  hY_p_nd->SetBarWidth(0.95);
  hY_p_nd->Draw("BAR");
  cout<< "Events in card 2 pedestal " <<EntriesCard2_p << endl;
  c1->cd(5);
  hY_p->GetXaxis()->SetTitle("Strip number");
  hY_p->GetYaxis()->SetTitleOffset(1.5);
  hY_p->GetYaxis()->SetTitle("Q (nC)");  
  hY_p->SetLineColor(2);
  hY_p->SetFillColor(2);
  hY_p->GetYaxis()->SetRangeUser(0,MaxYCont_p+MaxYCont_p*0.4);
  hY_p->SetBarWidth(0.95);
  hY_p->Draw("BAR");
  c1->cd(6);
  double *y_graph_p_y = new double[32];
  for(int i=0;i<32;i++) y_graph_p_y[i] =hY_p->GetBinContent(i+1);  
  TGraph* gY_p_mm = new TGraph(32, xaxis, y_graph_p_y);
  gY_p_mm->SetTitle("Strips Y - decoding applied");
  gY_p_mm->GetXaxis()->SetTitle("Detector side (mm)");
  gY_p_mm->GetYaxis()->SetTitle("Q (nC)");
  gY_p_mm->GetYaxis()->SetTitleOffset(1.5);
  gY_p_mm->GetYaxis()->SetRangeUser(0,MaxYCont_p+MaxYCont_p*0.4);
  gY_p_mm->SetLineColor(2);
  gY_p_mm->SetMarkerStyle(20);
  gY_p_mm->SetMarkerSize(1);
  gY_p_mm->Draw("APL");
  
 //===========================================
 //
 //      Treat the pedestal file  (comments as above)
 //                             
 //===========================================
  
  
      

 TFile*f_r = new TFile("../../5_12_2017/Run.root","read");  
 TTree*t_r = (TTree*)f_r->Get("DataTree");
 
 vector<double> *q_r           =0;
 vector<int>    *sat_r         =0;
 double          time_ns_r     =0;
 int             card_nb_r     =0;
   
 TBranch *charge_r     = 0;
 TBranch *saturation_r = 0;
 TBranch *card_r       = 0;
 TBranch *time_r       = 0;
 
 t_r->SetBranchAddress("time",&time_ns_r);
 t_r->SetBranchAddress("card",&card_nb_r);
 t_r->SetBranchAddress("Q",&q_r,&charge_r);
 t_r->SetBranchAddress("S",&sat_r,&saturation_r);
 
 
 TH1D*hX_r_nd = new TH1D("hX_r_nd","Strips X - decoding not applied", 32,0,32);                  // Create hitsograms of hX_p_nd = X strips, pedestal signal, not decoded (= decoding not applied)
 TH1D*hY_r_nd = new TH1D("hY_r_nd","Strips Y - decoding not applied", 32,0,32);                  // Create hitsograms of hY_p_nd = Y strips, pedestal signal, not decoded (= decoding not applied)

 TH1D*hX_r = new TH1D("hX_r","Strips X - decoding applied", 32,0,32);
 int EntriesCard1_r = 0;
 TH1D*hY_r = new TH1D("hY_r","Strips Y- decoding applied", 32,0,32);
 int EntriesCard2_r = 0;


 int nEntries_r=t_r->GetEntries(); 

 for (Int_t i = 0; i < nEntries_r; i++) {
      
      t_r->GetEntry(i);

      //cout << "time " << time_ns_r << endl;
      //cout << "card_nb " << card_nb_r << endl;


      Long64_t tentry_r = t_r->LoadTree(i);
      charge_r->GetEntry(tentry_r);
      saturation_r->GetEntry(tentry_r);
      //for (UInt_t j = 0; j < q->size(); ++j) {
      // printf("%lf   %d \n",  q->at(j) , sat->at(j));}
      if(card_nb_r==1)   {EntriesCard1_r++;
                          for (UInt_t j = 0; j < q_r->size(); ++j) {
                                              hX_r_nd->Fill(j, q_r->at(j));
                                              hX_r->Fill(strip_ydec[j]-1,q_r->at(j));}
                       }
 
      if(card_nb_r==2)   {EntriesCard2_r++;
                          for (UInt_t j = 0; j < q_r->size(); ++j) {
                                              hY_r_nd->Fill(j, q_r->at(j));
                                              hY_r->Fill(strip_ydec[j]-1,q_r->at(j));}
                       }

   }


   
   
   
   
   
     // Plots
  c2->cd(1);
  double MaxXCont_r = 0;
  for(int i = 0; i<32; i++) {if(hX_r_nd->GetBinContent(i)>MaxXCont_r) MaxXCont_r = hX_r_nd->GetBinContent(i);}
  hX_r_nd->GetXaxis()->SetTitle("Strip number");
  hX_r_nd->GetYaxis()->SetTitleOffset(1.5);
  hX_r_nd->GetYaxis()->SetTitle("Q (nC)");
  hX_r_nd->SetLineColor(2);
  hX_r_nd->SetFillColor(2);
  hX_r_nd->GetYaxis()->SetRangeUser(0,MaxXCont_r+MaxXCont_r*0.4);
  hX_r_nd->SetBarWidth(0.95);
  hX_r_nd->Draw("BAR");
  cout<<"Events in card 1 signal " << EntriesCard1_r << endl;
  c2->cd(2);
  hX_r->GetXaxis()->SetTitle("Strip number");
  hX_r->GetYaxis()->SetTitleOffset(1.5);
  hX_r->GetYaxis()->SetTitle("Q (nC)");
  hX_r->SetLineColor(2);
  hX_r->SetFillColor(2);
  hX_r->GetYaxis()->SetRangeUser(0,MaxXCont_r+MaxXCont_r*0.4);
  hX_r->SetBarWidth(0.95);
  hX_r->Draw("BAR");
  c2->cd(3);
  double *x_graph_r_y = new double[32];
  for(int i=0;i<32;i++) x_graph_r_y[i] =hX_r->GetBinContent(i+1);  
  TGraph* gX_r_mm = new TGraph(32, xaxis, x_graph_r_y);
  gX_r_mm->SetTitle("Strips X - decoding applied");
  gX_r_mm->GetXaxis()->SetTitle("Detector side (mm)");
  gX_r_mm->GetYaxis()->SetTitle("Q (nC)");
  gX_r_mm->GetYaxis()->SetTitleOffset(1.5);
  gX_r_mm->GetYaxis()->SetRangeUser(0,MaxXCont_r+MaxXCont_r*0.4);
  gX_r_mm->SetLineColor(2);
  gX_r_mm->SetMarkerStyle(20);
  gX_r_mm->SetMarkerSize(1);
  gX_r_mm->Draw("APL");
  
  
  
  c2->cd(4);
  double MaxYCont_r = 0;
  for(int i = 0; i<32; i++) {if(hY_r_nd->GetBinContent(i)>MaxYCont_r) MaxYCont_r= hY_r_nd->GetBinContent(i);}
  hY_r_nd->GetXaxis()->SetTitle("Strip number");
  hY_r_nd->GetYaxis()->SetTitleOffset(1.5);
  hY_r_nd->GetYaxis()->SetTitle("Q (nC)");
  hY_r_nd->SetLineColor(2);
  hY_r_nd->SetFillColor(2);
  hY_r_nd->GetYaxis()->SetRangeUser(0,MaxYCont_r+MaxYCont_r*0.4);
  hY_r_nd->SetBarWidth(0.95);
  hY_r_nd->Draw("BAR");
  cout<< "Events in card 2 signal " <<EntriesCard2_r << endl;
  c2->cd(5);
  hY_r->GetXaxis()->SetTitle("Strip number");
  hY_r->GetYaxis()->SetTitleOffset(1.5);
  hY_r->GetYaxis()->SetTitle("Q (nC)");  
  hY_r->SetLineColor(2);
  hY_r->SetFillColor(2);
  hY_r->GetYaxis()->SetRangeUser(0,MaxYCont_r+MaxYCont_r*0.4);
  hY_r->SetBarWidth(0.95);
  hY_r->Draw("BAR");
  c2->cd(6);
  double *y_graph_r_y = new double[32];
  for(int i=0;i<32;i++) y_graph_r_y[i] =hY_r->GetBinContent(i+1);  
  TGraph* gY_r_mm = new TGraph(32, xaxis, y_graph_r_y);
  gY_r_mm->SetTitle("Strips Y - decoding applied");
  gY_r_mm->GetXaxis()->SetTitle("Detector side (mm)");
  gY_r_mm->GetYaxis()->SetTitle("Q (nC)");
  gY_r_mm->GetYaxis()->SetTitleOffset(1.5);
  gY_r_mm->GetYaxis()->SetRangeUser(0,MaxYCont_r+MaxYCont_r*0.4);
  gY_r_mm->SetLineColor(2);
  gY_r_mm->SetMarkerStyle(20);
  gY_r_mm->SetMarkerSize(1);
  gY_r_mm->Draw("APL");

   
  
  //  Signal = Subtraction run/events - pedestal/events


 TH1D*hX_s_nd = new TH1D("hX_s","Strips X signal per event - decoding not applied", 32,0,32);
 TH1D*hY_s_nd = new TH1D("hY_s","Strips Y signal per event- decoding not applied", 32,0,32);

 TH1D*hX_s = new TH1D("hX_s","Strips X signal per event - decoding applied", 32,0,32);
 TH1D*hY_s = new TH1D("hY_s","Strips Y signal per event - decoding applied", 32,0,32);
 
 double *x_graph_s_y = new double[32];
 double *y_graph_s_y = new double[32];
 
 for (int i=0;i<32;i++) {hX_s_nd->Fill(i, (hX_r_nd->GetBinContent(i+1)/EntriesCard1_r)-(hX_p_nd->GetBinContent(i+1)/EntriesCard1_p));
                         hY_s_nd->Fill(i, (hY_r_nd->GetBinContent(i+1)/EntriesCard1_r)-(hY_p_nd->GetBinContent(i+1)/EntriesCard1_p));
                         hX_s->Fill(i, (hX_r->GetBinContent(i+1)/EntriesCard1_r)-(hX_p->GetBinContent(i+1)/EntriesCard1_p));
                         hY_s->Fill(i, (hY_r->GetBinContent(i+1)/EntriesCard1_r)-(hY_p->GetBinContent(i+1)/EntriesCard1_p));
                         x_graph_s_y[i] = (hX_r->GetBinContent(i+1)/EntriesCard1_r)-(hX_p->GetBinContent(i+1)/EntriesCard1_p);
                         y_graph_s_y[i] = (hY_r->GetBinContent(i+1)/EntriesCard1_r)-(hY_p->GetBinContent(i+1)/EntriesCard1_p);
                         }
                         
 c3->cd(1);
 hX_s_nd->GetXaxis()->SetTitle("Strip number");
 hX_s_nd->GetYaxis()->SetTitleOffset(1.5);
 hX_s_nd->GetYaxis()->SetTitle("Q (nC)");
 hX_s_nd->SetLineColor(2);
 hX_s_nd->SetFillColor(2);
 //hX_s_nd->GetYaxis()->SetRangeUser(0,MaxXCont_r+MaxXCont_r*0.4);
 hX_s_nd->SetBarWidth(0.95);
 hX_s_nd->Draw("BAR"); 
 
 c3->cd(2);
 hX_s->GetXaxis()->SetTitle("Strip number");
 hX_s->GetYaxis()->SetTitleOffset(1.5);
 hX_s->GetYaxis()->SetTitle("Q (nC)");
 hX_s->SetLineColor(2);
 hX_s->SetFillColor(2);
 //hX_s->GetYaxis()->SetRangeUser(0,MaxXCont_r+MaxXCont_r*0.4);
 hX_s->SetBarWidth(0.95);
 hX_s->Draw("BAR");

 c3->cd(3);
 TGraph* gX_s_mm = new TGraph(32, xaxis, x_graph_s_y);
 gX_s_mm->SetTitle("Strips X signal per event - decoding applied");
 gX_s_mm->GetXaxis()->SetTitle("Detector side (mm)");
 gX_s_mm->GetYaxis()->SetTitle("Q (nC)");
 gX_s_mm->GetYaxis()->SetTitleOffset(1.5);
 //gX_s_mm->GetYaxis()->SetRangeUser(0,MaxXCont_r+MaxXCont_r*0.4);
 gX_s_mm->SetLineColor(2);
 gX_s_mm->SetMarkerStyle(20);
 gX_s_mm->SetMarkerSize(1);
 gX_s_mm->Draw("APL");
 
 c3->cd(4);
 hY_s_nd->GetXaxis()->SetTitle("Strip number");
 hY_s_nd->GetYaxis()->SetTitleOffset(1.5);
 hY_s_nd->GetYaxis()->SetTitle("Q (nC)");
 hY_s_nd->SetLineColor(2);
 hY_s_nd->SetFillColor(2);
 //hX_s_nd->GetYaxis()->SetRangeUser(0,MaxXCont_r+MaxXCont_r*0.4);
 hY_s_nd->SetBarWidth(0.95);
 hY_s_nd->Draw("BAR"); 




 c3->cd(5);
 hY_s->GetXaxis()->SetTitle("Strip number");
 hY_s->GetYaxis()->SetTitleOffset(1.5);
 hY_s->GetYaxis()->SetTitle("Q (nC)");
 hY_s->SetLineColor(2);
 hY_s->SetFillColor(2);
 //hX_s_nd->GetYaxis()->SetRangeUser(0,MaxXCont_r+MaxXCont_r*0.4);
 hY_s->SetBarWidth(0.95);
 hY_s->Draw("BAR"); 
 
 c3->cd(6);
 TGraph* gY_s_mm = new TGraph(32, xaxis, y_graph_s_y);
 gY_s_mm->SetTitle("Strips Y signal per event - decoding applied");
 gY_s_mm->GetXaxis()->SetTitle("Detector side (mm)");
 gY_s_mm->GetYaxis()->SetTitle("Q (nC)");
 gY_s_mm->GetYaxis()->SetTitleOffset(1.5);
 //gX_s_mm->GetYaxis()->SetRangeUser(0,MaxXCont_r+MaxXCont_r*0.4);
 gY_s_mm->SetLineColor(2);
 gY_s_mm->SetMarkerStyle(20);
 gY_s_mm->SetMarkerSize(1);
 gY_s_mm->Draw("APL");

 

 //  f_r->Close();
  // f_p->Close();


  delete x_graph_p_y;
  delete y_graph_p_y;
  delete x_graph_r_y;
  delete y_graph_r_y;
  delete x_graph_s_y;
  delete y_graph_s_y;

}
*/