#define timewalk2bset2_cxx
#include "timewalk2bset2.h"
#include <TH2.h>
#include <TStyle.h>
#include <TCanvas.h>
#include <TH2.h>
#include <TStyle.h>
#include <TCanvas.h>
#include <TH2.h>
#include <TStyle.h>
#include <TCanvas.h>
#include <TH2.h>
#include <TStyle.h>
#include <TCanvas.h>
#include <TH2.h>
#include <TStyle.h>
#include <TCanvas.h>
#include <TH2.h>
#include <TStyle.h>
#include <TCanvas.h>
#include <TROOT.h>
#include <TChain.h>
#include <TFile.h>
#include <TH2.h>
#include <TStyle.h>
#include <TCanvas.h>
#include <TH2.h>
#include <TStyle.h>
#include <TCanvas.h>
#include <TH1.h>
#include <TH2.h>
#include <TStyle.h>
#include <TAttFill.h>
#include <TGraph.h>
#include <TFile.h>
#include <TAttMarker.h>
#include <TF1.h>
#include <TMath.h>
#include <TPad.h>
#include <TString.h>
#include <iostream>
#include <string>
#include <sstream>
#include <TCanvas.h>
#include <RooGlobalFunc.h>
#include <RooRealVar.h>
#include <RooGaussian.h>
#include <RooLandau.h>
#include <RooPlot.h>
#include <TCanvas.h>
#include <RooNumConvPdf.h>
#include <RooDataSet.h>
#include <RooDataHist.h>
#include  <TTree.h>
#include <TH1D.h>
#include <TRandom.h>
#include <RooGenericPdf.h>
#include <RooFFTConvPdf.h>
#include <TROOT.h>
#include <RooFitResult.h>
using namespace std;
using namespace RooFit;
//------------------------------------LANDAU FIT----------------------------------------------------
double fitlan(double *v,double *par)
   {
	   double x;
	    x=v[0];
	   double amp=par[0];
	   double mpv=par[1]+0.022278298;
	   double sigma=par[2];
	   
	   double result=amp*TMath::Landau(x,mpv,sigma);
	   return result;
   }

void timewalk2bset2::Loop(Int_t totevents)
{
//   In a ROOT session, you can do:
//      Root > .L timewalk2bset2.C
//      Root > timewalk2bset2 t
//      Root > t.GetEntry(12); // Fill t data members with entry number 12
//      Root > t.Show();       // Show values of entry 12
//      Root > t.Show(16);     // Read and show values of entry 16
//      Root > t.Loop();       // Loop on all entries
//

//     This is the loop skeleton where:
//    jentry is the global entry number in the chain
//    ientry is the entry number in the current Tree
//  Note that the argument to GetEntry must be:
//    jentry for TChain::GetEntry
//    ientry for TTree::GetEntry and TBranch::GetEntry
//
//       To read only selected branches, Insert statements like:
// METHOD1:
//    fChain->SetBranchStatus("*",0);  // disable all branches
//    fChain->SetBranchStatus("branchname",1);  // activate branchname
// METHOD2: replace line
//    fChain->GetEntry(jentry);       //read all branches
//by  b_branchname->GetEntry(ientry); //read only this branch
   if (fChain == 0) return;
   
    fio2b_set2=new TFile("timewalkcorrections2bset2.root","recreate");

   Long64_t nentries = fChain->GetEntriesFast();
    Long64_t nbytes = 0, nb = 0;
   // TH1F *h_risetime_1=new TH1F("risetime_10-90_adc1","risetime_1_check_adc1",5,0,3.5); //HISTOGRAM FOR FILLING RISETIME
    
    
   //for each event create a histogram
   for(int event=1;event<totevents+1;event++)
   {
	   hists_set2<<"histogram_adc1_set2"<<event;
	   h2b_set2_1=new TH1F(hists_set2.str().c_str(),"histogram_set2_1",101,-0.5,100.5);
	   hists_set2.str("");
	   
	   
	   int startEntry=1+100*(event-1);  //for event=1, the startEntry would be 1 so the for loop would start from 1. for the 2nd event the startEntry should be 101
	   int lastEntry=event*100;
	   
   // fill histogram with 100 entries	   
   for (Long64_t jentry=startEntry; jentry<lastEntry+1;jentry++) 
   {
      Long64_t ientry = LoadTree(jentry);
      if (ientry < 0) break;
      nb = fChain->GetEntry(jentry);   nbytes += nb;
      if(eventn_adc1>2000)
      h2b_set2_1->Fill(eventn_count,eventn_adc1_raw);
      
     // if (Cut(ientry) < 0) continue;
   }
       
	
       hists_offset<<"offset_set2"<<event;
	   hist_offset_set2=new TF1(hists_offset.str().c_str(),"pol0",0.5,100.5);
	   hists_offset.str("");
	   
	 
	    
	   fits_set2<<"histogram_set2_fits"<<event;
	   f_set2=new TF1(fits_set2.str().c_str(),fitlan,0.5,100.5,3);
	   fits_set2.str("");
	   
	  
      
		  
		  
			   cout<<"EVENT NUMBER= "<<event<<endl;
			    cout<<"**************************************************************************************"<<endl;
		       
		        
		         h2b_set2_1->Fit(hist_offset_set2,"Q","N",0.5,100.5);
                 h2b_set2_1->Add(hist_offset_set2,-1);
       
                 double binamplitude=h2b_set2_1->GetMaximumBin();
                 double x1=h2b_set2_1->GetXaxis()->GetBinCenter(binamplitude);
     
                 f_set2->SetParameter(0,binamplitude);
                 f_set2->SetParameter(1,x1);
                 f_set2->SetParameter(2,1.0);
    
                 h2b_set2_1->Fit(f_set2,"WL","N",40.5,70.5);
//--------------------------------------------------------RISETIME----------------------------------------------------------------------------
               /*   double y1_max=f_set2->GetMaximum(40.5,70.5);
                  // cout<<"the maximum amplitube value is y1= "<<y1_max<<endl;
                  double x=f_set2->GetX(y1_max,40.5,70.5);
                 // cout<<"the corresponding x-value is x= "<<x<<endl;
                  double y2=0.1*y1_max;
                // cout<<"the value of ten percent amplitude is y2= "<<y2<<endl;
                  double x2=f_set2->GetX(y2,40.5,53.5);
               // cout<<"the corresponding x value for 10% chek is x2= "<<x2<<endl;
                  double y3=0.9*y1_max;
              // cout<<"the value of 90% of the amplituvde is y3= "<<y3<<endl;
                  double x3=f_set2->GetX(y3,40.5,70.5);
              // cout<<"the corresponding x value for 90% check is x3= "<<x3<<endl;
     
                  double risetime_1=x3-x2;
             // cout<<"risetime check 10-90= "<<risetime_1<<endl;
                 double y4=0.3*y1_max;
            //cout<<"the value of 30% amplitude is y4= "<<y4<<endl;
                double x4=f_set2->GetX(y4,40.5,53.5);
           //cout<<"the corresponding x value for the 30% check is x4= "<<x4<<endl;
                double y5=0.7*y1_max;
            //cout<<"the value of 70% amplitude is y5= "<<y5<<endl;
               double x5=f_set2->GetX(y5,40.5,53.5);
           //cout<<"the corresponding x value is x5= "<<x5<<endl;
     
               double risetime_2=x5-x4;
               cout<<"risetime check 10-90 adc1= "<<risetime_1<<endl;
               cout<<"risetime check 30-70 adc1= "<<risetime_2<<endl;
     
    
                h_risetime_1->Fill(risetime_1); */
 
       
       
        
  // h_risetime_1->Draw("*L");     
   //write and delete the histogram  
   h2b_set2_1->Write();
   h2b_set2_1->Delete();
   fio2b_set2->Write();
}

}