Hi! I’m relatively new to using ROOT so I apologize if my question is naive. I’m attempting to fit two template histograms onto my data histogram with RooFit (I want to know how much of each template histogram RooFit thinks makes up my data histogram).
As far as I’ve learned by googling/reading the documentation, looks like the process is
histograms > RooDataHist’s > RooAddPdf (named model) and then model.fitTo(data)
I wrote all the code to do this, starting from two toyMC’s (h1 and h2), but I keep getting this error, along with my output stats box telling me it failed and that only about 15% of my data came from the 1st template (which is not true; as you can see in the code, it should be roughly 60%/40% split). This is the error I get (x100+):
RooAddPdf::model[ mc1_yield * modelmc1 + mc2_yield * modelmc2 ]
getLogVal() top-level p.d.f evaluates to zero @ !refCoefNorm=(), !pdfs=(modelmc1 = 0/9969,modelmc2 = 0/9986), !coefficients=(mc1_yield = 0.1,mc2_yield = 0.1)
Plotting the RooDataHist’s before the fit works fine, so I know the problem has to be with RooAddPdf or the fit. This is my code (if you want the toy MC templates I started with I can send or include those too, they are short):
[code]// This is a re-do of the toy MC example, but binned fit with RooFit
// Uses two exponential functions (‘templates’ h1 and h2) of 1mil data points each, from generate_out.root
// Generates roofit_toymc.root and associated histograms
// C Bertsche Spring 2013
// Note that RooFit tutorials can be found here: http://roofit.sourceforge.net/docs/tutorial/index.html
using namespace RooFit;
roofit_example2()
{
// S e t u p i n p u t f i l e s
// --------------------------------
// Start with original MC templates, make scaled-down MCs and ‘data’
TFile *setA = TFile::Open("generate_out.root");
TH1F *h1 = h1; // Get original MC templates (h1, h2)
TH1F *h2 = h2;
TH1F *h12 = new TH1F("h12"," ",500,0,10); // Create 'data' (h12)
TH1F *h1t = new TH1F("h1t"," ",500,0,10); // Create scaled MCs (h1t, h2t)
TH1F *h2t = new TH1F("h2t"," ",500,0,10);
for (int i=0;i<3000;i++) // Fill 'data' (h12) with one MC template
{
Double_t life=h1->GetRandom();
h12->Fill(life);
}
for (int i=0;i<2000;i++) // Fill 'data' (h12) with second MC template
{
Double_t life=h2->GetRandom();
h12->Fill(life);
}
for (int i=0;i<10000;i++) // Fill scaled MCs (h1t, h2t)
{
Double_t life1=h1->GetRandom();
h1t->Fill(life1);
Double_t life2=h2->GetRandom();
h2t->Fill(life2);
}
// Make RooFit histograms from input
RooRealVar x("x","x",0.0,10.0);
RooDataHist data("data","data",x,h12);
RooDataHist mc1("mc1","scaled mc1",x,h1t);
RooDataHist mc2("mc2","scaled mc2",x,h2t);
RooRealVar mc1_yield("mc1_yield","yield mc1",2000,0,5000); // These are variables for output
RooRealVar mc2_yield("mc2_yield","yield of mc2",2000,0,5000); // These are variables for output
// Make PDF from MC histograms
RooHistPdf modelmc1(“modelmc1”,“modelmc1”,x, mc1);
RooHistPdf modelmc2(“modelmc2”,“modelmc2”,x, mc2);
RooAddPdf model("model","model",RooArgList(modelmc1,modelmc2),RooArgList(mc1_yield, mc2_yield));
// Combines my MCs into one PDF model
// Plot the imported histogram(s)
RooPlot* dframe = x.frame(Title(“Data”));
data.plotOn(dframe);
RooPlot* mc1frame = x.frame(Title("MC Scaled (1)"));
mc1.plotOn(mc1frame);
RooPlot* mc2frame = x.frame(Title("MC Scaled (2)"));
mc2.plotOn(mc2frame);
TCanvas* c = new TCanvas("roofit_example","RooFit FractionFit Example: Input data/templates",800,1200);
c->Divide(1,3);
gROOT->SetStyle("Plain"); // Removes gray background from plots
c->cd(1) ; gPad->SetLeftMargin(0.15) ; dframe->GetYaxis()->SetTitleOffset(1.4) ; dframe->Draw();
c->cd(2) ; gPad->SetLeftMargin(0.15) ; mc1frame->GetYaxis()->SetTitleOffset(1.4) ; mc1frame->Draw();
c->cd(3) ; gPad->SetLeftMargin(0.15) ; mc2frame->GetYaxis()->SetTitleOffset(1.4) ; mc2frame->Draw();
// F i t a n d p l o t m o d e l
// -----------------------------------------------
// Now I have three RooDataHists: data, mc1, and mc2. They depend on RooRealVar x, x2, and x3 respectively
// The two MC RooDataHists have each been translated into a RooHistPdf (“modelmc1” and “modelmc2”)
// The two RooHistPdf’s have been combined into a RooAddPdf (“model”)
gROOT->SetStyle("Plain"); // removes gray background from plots
// RooAbsReal* nll = model.createNLL(data,NumCPU(1)) ; // From other example, looks like
// RooAbsReal* pll_phi = nll->createProfile(mc1_yield) ; // another way of doing the fitTo
// RooMinuit(*nll).migrad();
// RooMinuit(*nll).hesse();
// RooMinuit(*nll).minos(); //optional
model.fitTo(data, Extended ()); // It is this fitTo command that gives the statistical output
RooPlot* fitFrame=x.frame(Bins(50), Title("Fit Model"));
model.paramOn(fitFrame);
data.plotOn(fitFrame, RooFit::LineColor(kRed));
model.plotOn(fitFrame, LineStyle(kDashed));
model.plotOn(fitFrame, Components("modelmc1"), LineColor(kGreen));
model.plotOn(fitFrame, Components("modelmc2"), LineColor(kBlue));
TCanvas* c6 = new TCanvas("c6","Fit Model",800,1200);
gROOT->SetStyle("Plain"); // Removes gray background from plots
gPad->SetLeftMargin(0.15) ; fitFrame->GetYaxis()->SetTitleOffset(1.4) ; fitFrame->Draw();
RooArgSet* params = model.getParameters(x);
// O u t p u t r e s u l t s
// --------------------------
[/code]