#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 nRandom = 10000;
	int  model = 1;


	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", nRandom);
	hist_bkg -> FillRandom ("landau", nRandom);

	float scale_target_sig = 2.5;
	float scale_target_bkg = 1.1;
	hist_target -> FillRandom ("gaus",   nRandom*scale_target_sig);
	hist_target -> FillRandom ("landau", nRandom*scale_target_bkg);

	long nTarget = long (hist_target->Integral());
	hist_sig -> Scale(1/hist_sig->Integral());
	hist_bkg -> Scale(1/hist_bkg->Integral());


	// + Convert to RooFit
	//--------------------
	RooRealVar *varX    = new RooRealVar ("varX",    "x-variable", -5, 10);
	varX -> setBins (30);

	RooDataHist *roohist_target = new RooDataHist ("roohist_target", "conversion of target", *varX, RooFit::Import(*hist_target));
	RooDataHist *roohist_sig    = new RooDataHist ("roohist_sig",    "conversion of sig",    *varX, RooFit::Import(*hist_sig));
	RooDataHist *roohist_bkg    = new RooDataHist ("roohist_bkg",    "conversion of bkg",    *varX, RooFit::Import(*hist_bkg));

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

	RooRealVar *fracsig;
	RooRealVar *fracbkg;
	RooAddPdf  *fitmodel;

	if (model==1)
	{
		fracsig = new RooRealVar ("fracsig", "fraction of Sig", 0.5, 0.1, 0.9);
		fitmodel = new RooAddPdf (
				"fitmodel", "",
				RooArgList(*template_sig, *template_bkg),
				RooArgList(*fracsig));
	}
	else
	{
		fracsig = new RooRealVar ("fracsig", "fraction of Sig", 15000, 0.0, 30000);
		fracbkg = new RooRealVar ("fracbkg", "fraction of Bkg", 10000, 0.0, 30000);
		fitmodel = new RooAddPdf (
				"fitmodel", "",
				RooArgList(*template_sig, *template_bkg),
				RooArgList(*fracsig, *fracbkg));
	}

	fitmodel -> fitTo (*roohist_target);


	double  scale_sig;
	double  scale_bkg;
	if (model==1)
	{
		scale_sig = fracsig->getVal() * nTarget;
		scale_bkg = (1-fracsig->getVal()) * nTarget;
	}
	else
	{
		scale_sig = fracsig->getVal();
		scale_bkg = fracbkg->getVal();
	}

	printf (" Scale value for signal:     %.4f\n", scale_sig);
	printf (" Scale value for background: %.4f\n", scale_bkg);

	hist_sig -> Scale (scale_sig);
	hist_bkg -> Scale (scale_bkg);


	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(1);
	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");
}