RooFit Examples in CMSSW

Complete: 3

Introduction

This page is intended to provide examples with RooFit that are usable within CMSSW. The presented examples can be installed and compiled under CMSSW using the following recipe:

cmsrel CMSSW_3_1_0_pre7
cd CMSSW_3_1_0_pre7/src
cmsenv
addpkg PhysicsTools/CandExamples V02-00-00
cd PhysicsTools/CandExamples
scram b 
$CMSSW_BASE/test/slc4_ia32_gcc345/<name of the application<

Building RooFit applications in CMSSW

Before using the following examples, you should first understand how to build a RooFit application under CMSSW. This part is described in the following page:

Simple Maximum Likelihood fit with Toy Monte Carlo generation

Generate a Toy Monte Carlo sample according to a PDF modeled as Exponential plus a Gaussian, and perform an Extended Maximum Likelihood fit according to the same model.

Code fragments

Define the variable (x) and the parameters of interest (mu, sigma, lambda) with their validity range:
  RooRealVar x("x", "x", -10, 10);

  RooRealVar mu("mu", "average", 0, -1, 1);
  RooRealVar sigma("sigma", "sigma", 1, 0, 5);
  RooGaussian gauss("gauss","gaussian PDF", x, mu, sigma);
  RooRealVar lambda("lambda", "slope", -0.1, -5., 0.);

Define a Gaussian and an Exponential PDF according to the above parameters:

  RooGaussian gauss("gauss","gaussian PDF", x, mu, sigma);
  RooExponential expo("expo", "exponential PDF", x, lambda);

Define a linar combination of the two PDFs with yields s (signal) and b (background):

  RooRealVar s("s", "signal yield", 100, 0, 10000);
  RooRealVar b("b", "background yield", 100, 0, 10000);

  RooAddPdf sum("sum", "gaussian plus exponential PDF",
                RooArgList(gauss, expo), RooArgList(s, b));

Generate a mixed signal plus background sample with 2000 entries:

  RooDataSet * data = sum.generate(x, 2000);

Fit the sample according to the PDF model using Extended Unbinned Maximum Likelihood method:

  sum.fitTo(*data, RooFit::Extended());

Plot the fit result and save the output as EPS file:

  RooPlot * xFrame = x.frame() ;
  data->plotOn(xFrame) ;
  sum.plotOn(xFrame) ;
  sum.plotOn(xFrame, RooFit::Components(expo), RooFit::LineStyle(kDashed)) ;
  TCanvas c;
  xFrame->Draw();
  c.SaveAs("extendedLikelihoodFit.eps");

Result

Below the plot resulting from running the above example

extendedLikelihoodFit.gif

Complete example

Simple Chi-squared fit

Using the same PDF defined in the previous example, we can also adopt a binned chi-squared fit.

In order to defined the binned histogram, one has to use the RooDataHist class, setting the proper binning for the x variable:

  x.setBins(50);
  RooDataHist hist("hist", "hist", RooArgSet(x), *data);

Then, the chi-squared fit should be used instead of the likelihood fit. Notice, that if the yields should be fitted, an extra parameter (similar to the "extended" option in the likelihood fit) should be set to true in RooChi2Var:

  RooChi2Var chi2("chi2","chi2", sum, hist, true);
  RooMinuit minuit(chi2);
  minuit.migrad();
  minuit.hesse();

All the rest of the code remains unchanged.

Results

Below the plot resulting from running the above example

binnedChi2Fit.gif

Complete example

Multi-dimensional fits

RooFit does not distinguish, conceptually, multi-dimensional PDFs from one-dimensional ones. The only difference is that the support for plotting via RooPlot is limited to a single dimension, and in the case of multi-dimensional PDFs works with one-dimensional projectons.

The following example shows how to define a two-dimensional PDF using RooProdPdf as the product of independend one-dimensional PDFs, and how to project PDFs and data sets to a single dimension using RooPlot, or as two-dimensional standard ROOT contour plots:

  RooRealVar x("x", "x", -10, 10);
  RooRealVar y("y", "y", -10, 10);
  RooRealVar mux("mux", "average-x'", 0, -1, 1);
  RooRealVar sigmax("sigmax", "sigma-x'", 0.5, 0, 5);
  RooGaussian gaussx("gaussx","gaussian PDF x'", x, mux, sigmax);
  RooRealVar muy("muy", "average-y'", 0, -1, 1);
  RooRealVar sigmay("sigmay", "sigma-y'", 1.5, 0, 5);
  RooGaussian gaussy("gaussy","gaussian PDF y'", y, muy, sigmay);
  RooProdPdf gaussxy("gaussxy", "gaussxy", RooArgSet(gaussx, gaussy));

  RooRealVar lambdax("lambdax", "log-slope x", -0.1, -5., 0.);
  RooExponential expox("expox", "exponential PDF x", x, lambdax);
  RooRealVar lambday("lambday", "log-slope y", -0.1, -5., 0.);
  RooExponential expoy("expoy", "exponential PDF y", y, lambday);
  RooProdPdf expoxy("expoxy", "expoxy", RooArgSet(expox, expoy));

  RooRealVar s("s", "signal yield", 50, 0, 10000);
  RooRealVar b("b", "background yield", 200, 0, 10000);

  RooAddPdf sum("sum", "gaussian plus exponential PDF",
                RooArgList(gaussxy, expoxy), RooArgList(s, b));
  RooDataSet * data = sum.generate(RooArgSet(x,y), 250);
  RooFitResult* result = sum.fitTo(*data, RooFit::Extended(), RooFit::Minos(), RooFit::Save());
  assert(result != 0);
  result->Print();

  RooPlot * xFrame = x.frame() ;
  data->plotOn(xFrame) ;
  sum.plotOn(xFrame) ;
  sum.plotOn(xFrame, RooFit::Components(expoxy), RooFit::LineStyle(kDashed)) ;

  RooPlot * yFrame = y.frame() ;
  data->plotOn(yFrame) ;
  sum.plotOn(yFrame) ;
  sum.plotOn(yFrame, RooFit::Components(expoxy), RooFit::LineStyle(kDashed)) ;

  TCanvas c;
  xFrame->Draw();
  c.SaveAs("twoDimensionalFit-yProjection.eps");
  yFrame->Draw();
  c.SaveAs("twoDimensionalFit-xProjection.eps");
  TH2 * surf = dynamic_cast<TH2*>(sum.createHistogram("surf", x,
                                                      RooFit::Binning(100),
                                                      RooFit::YVar(y, RooFit::Binning(100))));
  surf->Draw("contz");
  c.SaveAs("twoDimensionalFit-xyContour.eps");
  TH2 * surfData = data->createHistogram(x, y, 100, 100, "", "surfData");
  surfData->Draw("contz");
  c.SaveAs("twoDimensionalFit-xyDataContour.eps");

Complete example

More examples

You can find more examples in the following pages:

Edit | Attach | Watch | Print version | History: r7 < r6 < r5 < r4 < r3 | Backlinks | Raw View | Raw edit | More topic actions...
Topic revision: r6 - 2009-05-27 - LucaLista
 
    • Cern Search Icon Cern Search
    • TWiki Search Icon TWiki Search
    • Google Search Icon Google Search

    CMSPublic All webs login

This site is powered by the TWiki collaboration platform Powered by PerlCopyright & 2008-2019 by the contributing authors. All material on this collaboration platform is the property of the contributing authors.
Ideas, requests, problems regarding TWiki? Send feedback