#include "TROOT.h"
#include "TRint.h"
#include "TH1F.h"
#include "TH2F.h"
#include "TF1.h"
#include "TCanvas.h"
#include "TChain.h"
#include "TVector3.h"
#include "TMath.h"
#include "TFile.h"
#include "TSystem.h"
#include "TGStatusBar.h"
#include "TSystem.h"
#include "TXMLEngine.h"
#include "TTree.h"
#include "TLorentzVector.h"
#include "TNtuple.h" 
#include <cstdlib>
#include<string> 
#include<iostream> 
#include<fstream>
#include<iomanip>
#include<cmath>
#include <vector>
#include "TStyle.h"
#include "TGraphErrors.h"
#include "TGraph.h"
#include "TLegend.h"
#include "TLine.h"

using namespace std;

int GetThetaBinned(Double_t);
int GetMomentaBinned(Double_t);
int GetZBinned(Double_t);
int GetPhiBinned(Double_t);

void BinnedPulls(){
   gStyle->SetOptFit(1011); 

   //-------------declare my histogram
   TH1D *Hist_[15][18][10][24];  //theta, mom, z, phi

   //-------------Bin numbers and range are defined here
   Double_t pt_lo = -3.0; Double_t pt_up = 3.0; Int_t pt_nbin = 200;
   //-------------strings to be pass to histo names
   TString PULL[1]={"PullOne_"};
   TString PHI[24]={"Phi_Neg180_165", "Phi_Neg165_150", "Phi_Neg150_135", "Phi_Neg135_120", "Phi_Neg120_105", "Phi_Neg105_90", "Phi_Neg90_75", "Phi_Neg75_60", "Phi_Neg60_45", "Phi_Neg45_30", "Phi_Neg30_15", "Phi_Neg15_0", "Phi_0_15", "Phi_15_30", "Phi_30_45", "Phi_45_60", "Phi_60_75", "Phi_75_90", "Phi_90_105", "Phi_105_120", "Phi_120_135", "Phi_135_150", "Phi_150_165", "Phi_165_180"};
   TString Z[10] = {"Z_Neg110", "Z_Neg110_105", "Z_Neg105_100", "Z_Neg100_95", "Z_Neg95_90", "Z_Neg90_85", "Z_Neg85_80", "Z_Neg80_75", "Z_Neg75_70", "Z_Neg70"};
   TString THETA[15]={"Theta0_10_","Theta10_15_","Theta15_20_","Theta20_25_","Theta25_30_","Theta35_40_","Theta40_45_","Theta45_50_","Theta50_60_","Theta60_70_","Theta70_80_","Theta80_90_","Theta90_110_","Theta110_","ThetaOther_"};
   TString MomentaBin[18]={"0_15", "15_30", "30_45", "45_60", "60_75", "75_90", "90_105", "105_120", "120_135", "135_150", "150_165", "165_180", "180_195", "195_210", "210_225", "225_240", "240_255", "255_ON"};
   //------------- Setting up Histgram Array

      for(int s=0;s<15;s++){
         for(int wj=0;wj<18;wj++){
            for (int mom =0; mom<10; mom++) {
               for (int cd =0; cd<24; cd++) {

                  char pt[180]; 
                  sprintf(pt,"%s%s%s%s%s%s","Hist_",THETA[s].Data(),MomentaBin[wj].Data(),Z[mom].Data(),PHI[cd].Data(),"\0"); 	  
                  Hist_[s][wj][mom][cd]=new TH1D(pt,pt,pt_nbin,pt_lo,pt_up);

 
               }//-Close Phi
            }//----Close Z
         }//------ Close mom
      }//----------Close theta

   
   TChain *fitter = new TChain("Fit");
   //fitter->Add("/Volumes/Mac_Storage/YIELD/SINGLE_TRACK/YIELD_DATA/PRODUCTION/RUN_56515/wcuts/small.root");
   fitter->Add("small.root");

   Int_t nEvent = fitter->GetEntries();  
   
   Double_t P_theta, Pull_One, Pull_Prob, P_p, P_Vz, P_phi;
   
   fitter->SetBranchAddress("P_theta",&P_theta);
   fitter->SetBranchAddress("P_phi",&P_phi);
   fitter->SetBranchAddress("P_Vz",&P_Vz);
   fitter->SetBranchAddress("Pull_One",&Pull_One);
   fitter->SetBranchAddress("P_p",&P_p);

   //------------Filling Histograms 
   for (Int_t i=0; i<nEvent; i++) {//nEvent
      if(!(i%500000)) std::cout << "done " << i << " out of " << nEvent << " ==> " << double(i)*100.0/double(nEvent) << "%" << std::endl;
      fitter->GetEntry(i);

      int PStat = GetThetaBinned(P_theta); 
      int Pmom_stat = GetMomentaBinned(P_p);     
      int PZ = GetZBinned(P_Vz);
      int PPhiStat = GetPhiBinned(P_phi); 
      if(Pull_Prob > 0.01){
         Hist_[PStat][Pmom_stat][PZ][PPhiStat]->Fill(Pull_One);
      }
   }
   
   TString PID[3]={"P_","Pip_","Pim_"};   
   TCanvas *c1 = new TCanvas("c1","Histo",800,650);
    
//################### HERE IS WHERE THE ANIMATED GIF IS CREATED ####################   
   for(int s=0;s<15;s++){    //--------------Theta
      cout<<"Making theta "<<s<<endl;
      
      char giffer[50];
      sprintf(giffer,"%s%s%s","ProtonMomPull_",THETA[s].Data(),".gif+200");
      //sprintf(giffer,"%s%s%s","/Volumes/Mac_Storage/YIELD/SINGLE_TRACK/YIELD_DATA/PRODUCTION/RUN_56515/FITKIN_Plots/P_/ProtonMomPull_",THETA[s].Data(),".gif+200");

      for(int wj=0;wj<18;wj++){  //-----------Momenta 18
         for (int mom =0; mom<10; mom++) { //--------------Z 10
            for (int cd =0; cd<24; cd++) { //--------------Phi
               c1->cd();
               Hist_[s][wj][mom][cd]->Draw();//Fit(g1,"R,Q");
               char axis[50];
               sprintf(axis,"%s%s%s",PID[0].Data(),"  ",PULL[0].Data()); 
               Hist_[s][wj][mom][cd]->GetXaxis()->SetTitle(axis);
               //sprintf(giffer,"%s%s%s%s%s%s","/Volumes/Mac_Storage/YIELD/SINGLE_TRACK/YIELD_DATA/PRODUCTION/RUN_56515/FITKIN_Plots/P_/ProtonMomPull_",THETA[s].Data(),MomentaBin[wj].Data(),Z[mom].Data(),PHI[cd].Data(),".jpeg");

               
               c1->Print(giffer);
               
            }//-Close Phi
         }//----Close Z
      }//------ Close mom
   }//----------Close theta
   
}
 

    
void smallHistoArrayII(){
  
   //Pulls();
   BinnedPulls();

    
}
  


int GetThetaBinned(Double_t theta){
  
  //Double_t theta = theta_in*180./3.14159;
  if(theta < 10.) return 0;
  else if(theta < 15) return 1;
  else if(theta < 20) return 2;
  else if(theta < 25) return 3;
  else if(theta < 30) return 4;
  else if(theta < 35) return 5;
  else if(theta < 40) return 6;
  else if(theta < 45) return 7;
  else if(theta < 50) return 8;
  else if(theta < 60) return 9;
  else if(theta < 70) return 10;
  else if(theta < 80) return 11;
  else if(theta < 90) return 12;
  else if(theta <110) return 13;
  else return 14;
}

int GetMomentaBinned(Double_t momenta){
   
   if(momenta < 0.15) return 0;
   else if(momenta < 0.3) return 1;
   else if(momenta < 0.45) return 2;
   else if(momenta < 0.60) return 3;
   else if(momenta < 0.75) return 4;
   else if(momenta < 0.9) return 5;
   else if(momenta < 1.05) return 6;
   else if(momenta < 1.2) return 7;
   else if(momenta < 1.35) return 8;
   else if(momenta < 1.5) return 9;
   else if(momenta < 1.65) return 10;
   else if(momenta < 1.8) return 11;
   else if(momenta < 1.95) return 12;
   else if(momenta <2.1) return 13;
   else if(momenta <2.25) return 14;
   else if(momenta <2.4) return 15;
   else if(momenta <2.55) return 16;   
   else return 17;
}

int GetZBinned(Double_t Zvert){
   if(Zvert < -110) return 0;
   else if(Zvert < -105) return 1;
   else if(Zvert < -100) return 2;
   else if(Zvert < -95) return 3;
   else if(Zvert < -90) return 4;
   else if(Zvert < -85) return 5;
   else if(Zvert < -80) return 6;
   else if(Zvert < -75) return 7;
   else if(Zvert < -70) return 8;   
   else return 9;
}

int GetPhiBinned(Double_t Phi){
   if(Phi < -165.) return 0;
   else if(Phi < -150.) return 1;
   else if(Phi < -135.) return 2;
   else if(Phi < -120.) return 3;
   else if(Phi < -105.) return 4;
   else if(Phi < -90.) return 5;
   else if(Phi < -75.) return 6;
   else if(Phi < -60.) return 7;
   else if(Phi < -45.) return 8;
   else if(Phi < -30.) return 9;
   else if(Phi < -15.) return 10;
   else if(Phi < 0.) return 11;
   else if(Phi < 15.) return 12;
   else if(Phi < 30) return 13;
   else if(Phi < 45) return 14;
   else if(Phi < 60) return 15;
   else if(Phi < 75) return 16;
   else if(Phi < 90) return 17;
   else if(Phi < 105) return 18;
   else if(Phi < 120) return 19;
   else if(Phi < 135) return 20;
   else if(Phi < 150) return 21;
   else if(Phi < 165) return 22;
   else return 23;


   
}
//   TString PID[3]={"P_","Pip_","Pim_"};
//   
//   TCanvas *c1 = new TCanvas("c1","Histo",800,650);
//   c1->Divide(2,1);
//   for(int ds=0; ds<2; ds++){//9
//      Int_t ID = (0+ds)/3;
//      char giffer[50];
//      sprintf(giffer,"%s%s%s%s%s%s","/Volumes/Mac_Storage/YIELD/SINGLE_TRACK/YIELD_DATA/PRODUCTION/RUN_56515/FITKIN_Plots/",PID[ID].Data(),"/Histogram_",PID[ID].Data(),PULL[ds].Data(),".gif+200"); 
//      for(int s=0;s<6;s++){ //6
//         for(int wj=0;wj<15;wj++){//15
//            for (int mom =0; mom <18; mom++) {//18
//               
//               c1->cd(1);
//               TF1 *g1 = new TF1("g1","gaus",pt_lo,pt_up); 
//               Hist_[ds][s][wj][mom]->Draw();//Fit(g1,"R,Q");
//               c1->cd(2);
//               Histcor_[ds][s][wj][mom]->Draw();//Fit(g1,"R,Q");
//               
//               char axis[50];
//               sprintf(axis,"%s%s%s",PID[ID].Data(),"  ",PULL[ds].Data()); 
//               Hist_[ds][s][wj][mom]->GetXaxis()->SetTitle(axis);
//               Histcor_[ds][s][wj][mom]->GetXaxis()->SetTitle(axis);
//               
////               char giffer[50];  //,"/",PULL[ds].Data()
////               sprintf(giffer,"%s%s%s%s%s%s%s%s%s","/Volumes/Mac_Storage/YIELD/SINGLE_TRACK/YIELD_DATA/PRODUCTION/RUN_56515/FITKIN_Plots/",PID[ID].Data(),"/Histogram_",PID[ID].Data(),PULL[ds].Data(),SEC[s].Data(),THETA[wj].Data(),MomentaBin[mom].Data(),".pdf");
//               
//               c1->Print(giffer);   
//               
//               
//            }//------Close Momenta
//         }//----------- Close Theta
//      }//-----------Close Sector
//   }//-----------Close Pull