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

#include "TH1D.h"
#include "RooRealVar.h"
#include "RooDataSet.h"
#include "RooGaussian.h"
#include "RooConstVar.h"
#include "RooAddPdf.h"
#include "RooChebychev.h"
#include "RooFitResult.h"
#include "TCanvas.h"
#include "TAxis.h"
#include "RooPlot.h"
#include "TFile.h"
#include "TStyle.h"
#include "TH2.h"


void fitter_RooFitQuestion ()
{
	long nTemplate = 10000;
	long nTarget   = 10000;
	// + Program starts
	//-----------------
	time_t time_01 = time(0);

	printf (" =================================================\n");
	printf ("  A program that fit 2 histograms to a target one\n");
	printf (" =================================================\n\n");

	TH1::AddDirectory(kFALSE);
	//gErrorIgnoreLevel = kError, kBreak, kSysError, kFatal;
	gStyle -> SetOptStat (0);



	// + Create histogram for fit test
	//--------------------------------
	// * For fit template
	int  nBin = 30;
	double bin_min = -5;
	double bin_max =  10;

	TH1D *hist_target = new TH1D ("hist_target", "", nBin, bin_min, bin_max);
	TH1D *hist_sig    = new TH1D ("hist_sig",    "", nBin, bin_min, bin_max);
	TH1D *hist_bkg    = new TH1D ("hist_bkg",    "", nBin, bin_min, bin_max);
	hist_target -> Sumw2();
	hist_sig    -> Sumw2();
	hist_bkg    -> Sumw2();


	// + Generate event for template
	//------------------------------
	hist_sig -> FillRandom ("gaus", nTemplate/2);
	hist_bkg -> FillRandom ("landau", nTemplate/5);
	hist_target -> FillRandom ("gaus", nTarget);
	hist_target -> FillRandom ("landau", nTarget/3);


	// + Convert to RooFit
	//--------------------
	RooRealVar *varX    = new RooRealVar ("varX",    "x-variable", -5, 5);
	RooRealVar *fracsig = new RooRealVar ("fracsig", "fraction of signal", 0.3, 0.0, 1.0);

	RooDataHist *roohist_target = new RooDataHist ("roohist_target", "conversion of hist_target", *varX, hist_target);
	RooDataHist *roohist_sig    = new RooDataHist ("roohist_sig",    "conversion of hist_sig",    *varX, hist_sig);
	RooDataHist *roohist_bkg    = new RooDataHist ("roohist_bkg",    "conversion of hist_bkg",    *varX, hist_bkg);

	RooHistPdf *template_sig = new RooHistPdf ("template_sig", "histogram-pdf for signal",     RooArgSet(*varX), *roohist_sig);
	RooHistPdf *template_bkg = new RooHistPdf ("template_bkg", "histogram-pdf for background", RooArgSet(*varX), *roohist_bkg);

	RooAddPdf *fitmodel = new RooAddPdf ("fitmodel", "", RooArgList(*template_sig, *template_bkg), RooArgList(*fracsig));

	//RooFitResult *fitresult = (RooFitResult*)fitmodel -> fitTo (*roohist_target, Save());
	RooFitResult *fitresult = (RooFitResult*)fitmodel -> fitTo (*roohist_target);

	RooRealVar *roofrac_sig = (RooRealVar*)fitresult -> floatParsFinal().find("fracsig");

	double  fitvalue_sigfrac = roofrac_sig -> getVal();
	double  scale_sig = fitvalue_sigfrac * hist_target->Integral();
	double  scale_bkg = (1-fitvalue_sigfrac) * hist_target->Integral();

	hist_sig -> Scale (scale_sig);
	hist_bkg -> Scale (scale_bkg);
	//hist_sig -> Scale (2);
	//hist_bkg -> Scale (5.0/3.0);


	hist_target -> SetMarkerStyle (20);
	hist_target -> SetMarkerSize  (0.8);
	hist_target -> SetMarkerColor (kBlack);

	hist_sig -> SetFillColor (kOrange);
	hist_bkg -> SetFillColor (kAzure+2);

	THStack *stack_afterfit = new THStack ("stack_afterfit", "");
	stack_afterfit -> Add (hist_bkg);
	stack_afterfit -> Add (hist_sig);

	TCanvas *canvas = new TCanvas ("canvas", "", 800, 800);
	canvas -> cd();

	hist_target -> Draw ("ep");
	stack_afterfit -> Draw ("hist same");
	hist_target -> Draw ("ep same");


	system ("mkdir -p  Output/PlotRooFit/");
	canvas -> SaveAs ("Output/PlotRooFit/test.png");
}