#include <inttypes.h>
#include <iostream>
#include <fstream>
#include <utility>
#include <vector>
#include <map>
#include <stdio.h>
#include <stdlib.h>
#include <string>
#include <cstring>
#include <sstream>

#include "TLatex.h"
#include "TPaveStats.h"
#include "TGraphPainter.h"
#include "TString.h"
#include "TCollection.h"
#include "TCanvas.h"
#include "TFile.h"
#include "TH1F.h"
#include "TF1.h"
#include "TGraph.h"
#include "TGraphErrors.h"
#include "TMinuit.h"
#include "TKey.h"
#include "TDatime.h"
#include "TAxis.h"
#include "TLine.h"
#include "TTree.h"
#include "TBranch.h"
#include "TStyle.h"
#include "TROOT.h"
#include "TSystem.h"

Double_t propo_func (Double_t *x, Double_t *par) {
    
    // Returns funcition y = Ax
    
    Double_t slope = par[0];            // slope
    Double_t x1 = x[0] ;            // input pulse height in counts
    
    // cout << " x1=" << x1 << " slope=" << slope << endl;
    
    if (x1 <= 0) return 0;
    
    return slope*x1;
}
void SiPM_waveform_integral(void) {

	//gStyle->Reset();
	//gROOT->SetStyle("Plain");
    gStyle->SetPalette(1,0);
    // gStyle->SetOptStat(kTRUE);
    gStyle->SetOptStat(kFALSE);
    // gStyle->SetOptFit(kTRUE);
    //gStyle->SetOptFit(kFALSE);
    // gStyle->SetOptFit(11111);
    // gStyle->SetOptStat(2220);
	/*gStyle->SetStatH(0.10);
	gStyle->SetStatW(0.60);
	gStyle->SetStatX(0.99);
	gStyle->SetStatY(0.99);*/

    gStyle->SetPadRightMargin(0.15);
    gStyle->SetPadLeftMargin(0.2);
    gStyle->SetPadBottomMargin(0.15);
    //gStyle->SetFillColor(0);
    gStyle->SetTitleOffset(1.2, "Y");
  	gStyle->SetTitleSize(0.07,"xyz");
	gStyle->SetTitleSize(0.08,"h");
	gStyle->SetLabelSize(0.07,"xyz");
	gStyle->SetTitleX(0);
	gStyle->SetTitleAlign(13);
	gStyle->SetNdivisions(505,"xy");
    
    Double_t xmin, xmax, ymin, ymax;
    Double_t slope=1;
    TF1 *propo = new TF1 ("propo",propo_func,xmin,xmax,1);
    propo->SetParameter(0,slope);
    propo->SetParName(0,"slope");
    
    // get histograms
    
    TString filename = "waveform_dump_R566_n100.root";
    TFile *file = new TFile(filename,"read");
    if (file != NULL) {
        cout << " File opened=" << filename.Data() << endl;
    }
    else {
        cout << " Could not open file=" << filename.Data() << endl;
        return;
    }
    
    const Int_t nhist=100;
    TH1F *wavehist[nhist];
    for (Int_t kk=0; kk<nhist; kk++) {
        TString hname;
        hname.Form("waveform%03d",kk+3);               // names begin from event 3
        wavehist[kk] = (TH1F*)file->Get(hname.Data());
        // cout << " Get histogram=" << hname.Data() << endl;
    }
    
    // create TProfile histogram of integrals
    
    const Int_t nsamples=wavehist[0]->GetNbinsX();
    TProfile *p_int = new TProfile("p_int","Integral vs sample",nsamples,0,nsamples);
    
    for (Int_t kk=0; kk<nhist; kk++) {
        // pedestal
        Double_t pedestal = 0;
        cout << " Process histogram = " << kk << endl;
        for (Int_t jj=0; jj<10; jj++) {
            pedestal += wavehist[kk]->GetBinContent(jj+1);
            // cout << " jj=" << jj << " kk=" << kk << " pedestal=" << pedestal << endl;
        }
        pedestal /= 10;
        cout << " kk=" << kk << " pedestal=" << pedestal << endl << endl;
        
        // contents
        Double_t integ = 0;
        for (Int_t jj=10; jj<nsamples; jj++) {
            Double_t content = wavehist[kk]->GetBinContent(jj+1) - pedestal;
            // cout << " jj=" << jj << " kk=" << kk << " content=" << content << endl;
            integ += content;
            p_int->Fill(jj,integ);
        }
        cout << endl;
    }
    
    // scale to maximum
    
    p_int->Scale(1/p_int->GetMaximum());
    
    TCanvas *c0 = new TCanvas("c0","c0 SiPM_waveform_integral",200,10,1000,500);
    c0->Divide(2,1);
    
    c0->cd(1);
    gPad->SetGridx();
    gPad->SetGridy();
    
    
    ymin=0;
    ymax=2;
    xmin=0;
    xmax=100;
    wavehist[0]->SetTitle("");
    wavehist[0]->GetYaxis()->SetLabelSize(0.05);
    wavehist[0]->GetYaxis()->SetTitleSize(0.07);
    wavehist[0]->GetXaxis()->SetLabelSize(0.05);
    wavehist[0]->GetXaxis()->SetTitleSize(0.07);
    wavehist[0]->GetYaxis()->SetTitleOffset(1.5);
    wavehist[0]->GetYaxis()->SetTitle("FADC (counts)");
    wavehist[0]->GetXaxis()->SetTitle("FADC sample (4 ns)");
    wavehist[0]->GetXaxis()->SetNdivisions(505);
    wavehist[0]->SetMinimum(0);
    wavehist[0]->SetLineColor(2);
    // wavehist[0]->GetXaxis()->SetRangeUser(xmin,xmax);
    // wavehist[0]->GetYaxis()->SetRangeUser(ymin,ymax);
    wavehist[0]->SetMarkerStyle(20);
    wavehist[0]->SetMarkerColor(2);
    wavehist[0]->Draw();
    
    TLine *line = new TLine(0,100,100,100);
    line->Draw();
    
    
    // TCanvas *c1 = new TCanvas("c1","c1 SiPM_waveform_integral",200,10,800,700);
    
    c0->cd(2);
    gPad->SetGridx();
    gPad->SetGridy();
    
    ymin=0;
    ymax=1.5;
    xmin=0;
    xmax=100;
    p_int->SetTitle("");
    p_int->GetYaxis()->SetLabelSize(0.05);
    p_int->GetYaxis()->SetTitleSize(0.07);
    p_int->GetXaxis()->SetLabelSize(0.05);
    p_int->GetXaxis()->SetTitleSize(0.07);
    p_int->GetYaxis()->SetTitleOffset(1.5);
    p_int->GetYaxis()->SetTitle("Normalized Integral");
    p_int->GetXaxis()->SetTitle("FADC sample (4 ns)");
    p_int->GetXaxis()->SetNdivisions(505);
    p_int->SetMinimum(0);
    p_int->SetLineColor(2);
    // p_int->GetXaxis()->SetRangeUser(xmin,xmax);
    // p_int->GetYaxis()->SetRangeUser(ymin,ymax);
    p_int->SetMarkerStyle(20);
    p_int->SetMarkerColor(2);
    p_int->SetMarkerSize(0.4);
    p_int->Draw();

    
    
    /*TLegend *leg = new TLegend(0.2,0.8,0.5,0.9);
    leg->AddEntry(gr_SiPM_int,"Data","p");
    leg->AddEntry(gr_SiPM_corrected,"Data Corrected","p");
    leg->Draw();
    
    TString text;
    text.Form("Correction factors");
    TLatex *t3 = new TLatex(0.25,0.7,text);
    t3->SetNDC();
    t3->SetTextSize(0.03);
    t3->Draw();
    text.Form("k=%.3f",k);
    t3->DrawLatex(0.25,0.66,text);
    text.Form("Npixels=%.0f",Npixels);
    t3->DrawLatex(0.25,0.62,text);
    
    text.Form("Flux scale=%.3f",flux_fudge);
    t3->DrawLatex(0.25,0.58,text);*/
    
    c0->SaveAs("SiPM_waveform_integral.png");
    c0->SaveAs("SiPM_waveform_integral.pdf");
}