#include "TFile.h"
#include "TCanvas.h"
#include "TStyle.h"
#include "TH1.h"
#include "TH2.h"
#include "TH3.h"
#include "TGaxis.h"
#include "TRandom.h"
#include "TLegend.h"
#include "TPaveStats.h"
#include "TGraph.h"
#include "TPad.h"

float offx=1.2;
float offy=1.0;
double Emin=1200;
double Emax=5000;
int BinSimMin0=0;
int BinSimMax0=0;
int BinSimMin097=0;
int BinSimMax097=0;
int BinSimMin0sum=0;
int BinSimMax0sum=0;
int BinExpMin0=0;
int BinExpMax0=0;
int BinBibMin0=0;
int BinBibMax0=0;
int BinBibMin1=0;
int BinBibMax1=0;
int BinSimMin1=0;
int BinSimMax1=0;
int BinSimMin197=0;
int BinSimMax197=0;
int BinSimMin1sum=0;
int BinSimMax1sum=0;
int BinExpMin1=0;
int BinExpMax1=0;
float margr=0.08;
   float w=3;
   float margl=0.12;

void detsumbibratio()
{
	gStyle->SetOptStat(00000);
	gStyle->SetOptFit(0);
		TFile *f = TFile::Open("SimRes.root");
   if (f == 0) {
      // if we cannot open the file, print an error message and return immediatly
      printf("Error: cannot open SimLuna 6998keV.root!\n");
      return;
   }
   f->ls();
   TFile *g = TFile::Open("run0611_GePD.root");
   if (g == 0) {
      // if we cannot open the file, print an error message and return immediatly
      printf("Error: cannot open Experimental GePD File!\n");
      return;
   }
    g->ls();
     TFile *h = TFile::Open("run0611_GeDD.root");
   if (h == 0) {
      // if we cannot open the file, print an error message and return immediatly
      printf("Error: cannot open Experimental GeDD File!\n");
      return;
   }
    h->ls();
    TFile *i = TFile::Open("SimRes.root");
   if (i == 0) {
      // if we cannot open the file, print an error message and return immediatly
      printf("Error: cannot open SimLuna 6997keV.root!\n");
      return;
   }
   i->ls();
   TFile *l = TFile::Open("run0735.root");
   if (l == 0) {
      // if we cannot open the file, print an error message and return immediatly
      printf("Error: cannot open Bib.root!\n");
      return;
   }
   l->ls();
   	TTree *t0=0;
   	TTree *t097=0;
   		f->GetObject("Tree1",t0);
			if (t0 == nullptr ) {
    			  printf("Error: cannot get Tree1!\n");
      			return;
   			}
   			i->GetObject("Tree1",t097);
			if (t097 == nullptr ) {
    			  printf("Error: cannot get Tree1!\n");
      			return;
   			}
    TCanvas *c01 = new TCanvas("c01","c01",1280,1024);
	TPad *pad0a = new TPad("pad0a", "pad0a", 0.0, 0.3, 1.0, 1.0);
	pad0a->Draw();             // Draw the upper pad: pad0a
   pad0a->cd();
   pad0a->SetLogy();
    gPad->SetLeftMargin(margl);
   gPad->SetRightMargin(margr);
	gROOT->cd(); // newly created histograms should go here
	 TH1F *hExp0 = (TH1F *)g->Get("GePD");
	 TH1F *hBib0 = (TH1F *)l->Get("GePD");
        TH1F *hSim0 = (TH1F*)hExp0->Clone("htemp0");
        hSim0->Reset("M");
		TString hString0 = TString::Format("Edep[0] >> htemp0"); 
			t0->Draw(hString0);	
   			BinSimMin0 = hSim0->GetXaxis()-> FindBin(Emin);
   			BinSimMax0 = hSim0->GetXaxis()-> FindBin(Emax);
   			hSim0->Scale(1./hSim0->Integral(BinSimMin0,BinSimMax0), "width");
   		hSim0->GetXaxis()->SetRangeUser(0,7500);
   			hSim0->SetOption("HIST"); 
			c01->Modified(); 
			c01->Update(); 
   			hSim0->SetTitle("; Energy (keV); Counts");
   			hSim0->SetFillColor(kBlue);
   			hSim0->SetLineColor(kBlue);
   			hSim0->SetFillStyle(3005);
   			c01->Update();
   			TH1F *hSim097 = (TH1F*)hExp0->Clone("htemp097");
        hSim097->Reset("M");
		TString hString097 = TString::Format("Edep[0] >> htemp097"); 
			t097->Draw(hString097,"","same");	
   			BinSimMin097 = hSim097->GetXaxis()-> FindBin(Emin);
   			BinSimMax097 = hSim097->GetXaxis()-> FindBin(Emax);
   			hSim097->Scale(1./hSim097->Integral(BinSimMin097,BinSimMax097), "width");
   		  hSim097->GetXaxis()->SetRangeUser(0,7500);
   			hSim097->SetOption("HIST"); 
			c01->Modified(); 
			c01->Update(); 
   			hSim097->SetTitle("; Energy (keV); Counts");
   			hSim097->SetFillColor(kMagenta);
   			hSim097->SetLineColor(kMagenta);
   			hSim097->SetFillStyle(3005);
   			c01->Update();
   			TH1F *hSim0sum = (TH1F *)hSim0->Clone("hSim0sum");
    		if (!(hSim0sum->GetSumw2N() > 0))hSim0sum->Sumw2(kTRUE);
			hSim0sum->Add(hSim097,1);
			hSim0sum->GetXaxis()->SetRangeUser(0,7500);
   			hSim0sum->SetOption("HIST"); 
   			c01->Update();
   			hSim0sum->Draw("HIST SAME");
        c01->Update();
			c01->Modified(); 
			c01->Update(); 
   			hSim0sum->SetTitle("; Energy (keV); Counts");
   			hSim0sum->SetFillColor(kOrange);
   			hSim0sum->SetLineColor(kOrange);
   			hSim0sum->SetFillStyle(3005);
   			c01->Update();
   				BinExpMin0 = hExp0->GetXaxis()-> FindBin(Emin);
   			BinExpMax0 = hExp0->GetXaxis()-> FindBin(Emax);
        hExp0->Scale(1./hExp0->Integral(BinExpMin0,BinExpMax0), "width");
       hExp0->SetOption("HIST"); 
		c01->Modified(); 
		c01->Update(); 
        hExp0->SetLineColor(kRed);
        hExp0->SetFillColor(kRed);
        hExp0->SetFillStyle(3005);
        hExp0->SetTitle("; Energy (keV); Counts");
        hExp0->GetXaxis()->SetTitleOffset(offx);
        hExp0->GetYaxis()->SetTitleOffset(offy);
        hExp0->SetStats(0);
   		c01->Update();
        hExp0->Draw("HIST SAME");
        c01->Update();
        	BinBibMin0 = hBib0->GetXaxis()-> FindBin(Emin);
   			BinBibMax0 = hBib0->GetXaxis()-> FindBin(Emax);
        hBib0->Scale(1./hBib0->Integral(BinBibMin0,BinBibMax0), "width");
       hBib0->SetOption("HIST"); 
		c01->Modified(); 
		c01->Update(); 
        hBib0->SetLineColor(kBlack);
        hBib0->SetFillColor(kBlack);
        hBib0->SetFillStyle(3005);
        hBib0->SetTitle("; Energy (keV); Counts");
        hBib0->GetXaxis()->SetTitleOffset(offx);
        hBib0->GetYaxis()->SetTitleOffset(offy);
        hBib0->SetStats(0);
   		c01->Update();
        hBib0->Draw("HIST SAME");
        c01->Update();
       	TH1F *hSimBib0sum = (TH1F *)hSim0sum->Clone("hSimBib0sum");
    		if (!(hSimBib0sum->GetSumw2N() > 0))hSimBib0sum->Sumw2(kTRUE);
			hSimBib0sum->Add(hBib0,1);
			hSimBib0sum->GetXaxis()->SetRangeUser(0,7500);
   			hSimBib0sum->SetOption("HIST"); 
   			c01->Update();
   			hSimBib0sum->Draw("HIST SAME");
        c01->Update();
			c01->Modified(); 
			c01->Update(); 
   			hSimBib0sum->SetTitle("; Energy (keV); Counts");
   			hSimBib0sum->SetFillColor(kViolet);
   			hSimBib0sum->SetLineColor(kViolet);
   			hSimBib0sum->SetFillStyle(3005);
   			c01->Update();
        	TLegend* leg0 = new TLegend(0.35, 0.7, .5, .8);
   			leg0->SetHeader("");
  			leg0->SetNColumns(1);
   			leg0->AddEntry(hSim0, "MC 6998 keV", "l");
   			leg0->AddEntry(hSim097, "MC 6997 keV", "l");
   			leg0->AddEntry(hSim0sum, "MC", "l");
   			leg0->AddEntry(hExp0, "Data: Run0611", "l");
   			leg0->AddEntry(hBib0, "Data BIB: Run0735", "l");
   			leg0->AddEntry(hSimBib0sum, "MC+BIB", "l");
   			 leg0->SetBorderSize(0);
   			leg0->Draw(); //if wanted
       
      c01->cd();
   c01->Update();
   TPad *pad0b = new TPad("pad0b", "pad0b", 0.0, 0.0, 1.0, 0.3);
    pad0b->SetTopMargin(0.0);
   pad0b->SetGridx(); // vertical grid
   pad0b->SetLogy();
   pad0b->Draw();
   pad0b->cd();       // pad0b becomes the current pad
    TH1F *ratio0 = (TH1F *)hExp0->Clone("ratio0");
    if (!(ratio0->GetSumw2N() > 0))ratio0->Sumw2(kTRUE);
	ratio0->Divide(hSim0sum);
	ratio0-> SetLineColor(kGreen);
	ratio0->GetXaxis()->SetRangeUser(0,7500);
	ratio0->SetFillColor(kGreen);
	ratio0->SetTitle("; Energy (keV); ratio");
	ratio0->Draw("EP");
	TLegend* leg0r = new TLegend(0.45, 0.8, .6, .9);
   	leg0r->SetHeader("");
    leg0r->SetNColumns(1);
   	leg0r->AddEntry(ratio0, "Data/MC", "l");
   	leg0r->SetBorderSize(0);
   	leg0r->Draw(); //if wanted

   gPad->SetLeftMargin(margl);
   gPad->SetRightMargin(margr);
   gPad->SetBottomMargin(0.15);
   
       c01->Print("DepEne_SimSum69976998ExpBibRatioScale12005000_Ep=278_04__GePD.png"); 
      delete c01;
   
}