Organization and simultaneous fits: easy interactive access to workspace contents - CINT to CLING code migration
Author: Wouter Verkerke
This notebook tutorial was automatically generated with ROOTBOOK-izer from the macro found in the ROOT repository on Tuesday, March 19, 2024 at 07:16 PM.
using namespace RooFit;
void fillWorkspace(RooWorkspace &w);
Definition of a helper function:
%%cpp -d
void fillWorkspace(RooWorkspace &w)
{
// C r e a t e p d f a n d f i l l w o r k s p a c e
// --------------------------------------------------------
// Declare observable x
RooRealVar x("x", "x", 0, 10);
// Create two Gaussian PDFs g1(x,mean1,sigma) anf g2(x,mean2,sigma) and their parameters
RooRealVar mean("mean", "mean of gaussians", 5, 0, 10);
RooRealVar sigma1("sigma1", "width of gaussians", 0.5);
RooRealVar sigma2("sigma2", "width of gaussians", 1);
RooGaussian sig1("sig1", "Signal component 1", x, mean, sigma1);
RooGaussian sig2("sig2", "Signal component 2", x, mean, sigma2);
// Build Chebychev polynomial pdf
RooRealVar a0("a0", "a0", 0.5, 0., 1.);
RooRealVar a1("a1", "a1", 0.2, 0., 1.);
RooChebychev bkg("bkg", "Background", x, RooArgSet(a0, a1));
// Sum the signal components into a composite signal pdf
RooRealVar sig1frac("sig1frac", "fraction of component 1 in signal", 0.8, 0., 1.);
RooAddPdf sig("sig", "Signal", RooArgList(sig1, sig2), sig1frac);
// Sum the composite signal and background
RooRealVar bkgfrac("bkgfrac", "fraction of background", 0.5, 0., 1.);
RooAddPdf model("model", "g1+g2+a", RooArgList(bkg, sig), bkgfrac);
w.import(model);
}
Create a workspace named 'w' With CINT w could exports its contents to a same-name C++ namespace in CINT 'namespace w'. but this does not work anymore in CLING. so this tutorial is an example on how to change the code
RooWorkspace *w1 = new RooWorkspace("w", true);
Fill workspace with pdf and data in a separate function
fillWorkspace(*w1);
[#0] WARNING:InputArguments -- The parameter 'sigma1' with range [-inf, inf] of the RooGaussian 'sig1' exceeds the safe range of (0, inf). Advise to limit its range. [#0] WARNING:InputArguments -- The parameter 'sigma2' with range [-inf, inf] of the RooGaussian 'sig2' exceeds the safe range of (0, inf). Advise to limit its range. [#1] INFO:ObjectHandling -- RooWorkspace::import(w) importing RooAddPdf::model [#1] INFO:ObjectHandling -- RooWorkspace::import(w) importing RooChebychev::bkg [#1] INFO:ObjectHandling -- RooWorkspace::import(w) importing RooRealVar::x [#1] INFO:ObjectHandling -- RooWorkspace::import(w) importing RooRealVar::a0 [#1] INFO:ObjectHandling -- RooWorkspace::import(w) importing RooRealVar::a1 [#1] INFO:ObjectHandling -- RooWorkspace::import(w) importing RooRealVar::bkgfrac [#1] INFO:ObjectHandling -- RooWorkspace::import(w) importing RooAddPdf::sig [#1] INFO:ObjectHandling -- RooWorkspace::import(w) importing RooGaussian::sig1 [#1] INFO:ObjectHandling -- RooWorkspace::import(w) importing RooRealVar::mean [#1] INFO:ObjectHandling -- RooWorkspace::import(w) importing RooRealVar::sigma1 [#1] INFO:ObjectHandling -- RooWorkspace::import(w) importing RooRealVar::sig1frac [#1] INFO:ObjectHandling -- RooWorkspace::import(w) importing RooGaussian::sig2 [#1] INFO:ObjectHandling -- RooWorkspace::import(w) importing RooRealVar::sigma2
Print workspace contents
w1->Print();
RooWorkspace(w) w contents variables --------- (a0,a1,bkgfrac,mean,sig1frac,sigma1,sigma2,x) p.d.f.s ------- RooChebychev::bkg[ x=x coefList=(a0,a1) ] = 0.8 RooAddPdf::model[ bkgfrac * bkg + [%] * sig ] = 0.9/1 RooAddPdf::sig[ sig1frac * sig1 + [%] * sig2 ] = 1/1 RooGaussian::sig1[ x=x mean=mean sigma=sigma1 ] = 1 RooGaussian::sig2[ x=x mean=mean sigma=sigma2 ] = 1
this does not work anymore with CLING use normal workspace functionality
use normal workspace methods
RooAbsPdf *model = w1->pdf("model");
RooRealVar *x = w1->var("x");
std::unique_ptr<RooDataSet> d{ model->generate(*x, 1000)};
std::unique_ptr<RooFitResult> r{model->fitTo(*d, PrintLevel(-1))};
[#1] INFO:Fitting -- RooAbsPdf::fitTo(model) fixing normalization set for coefficient determination to observables in data [#1] INFO:Fitting -- using CPU computation library compiled with -mavx2 [#1] INFO:Fitting -- RooAddition::defaultErrorLevel(nll_model_modelData) Summation contains a RooNLLVar, using its error level [#1] INFO:Minimization -- RooAbsMinimizerFcn::setOptimizeConst: activating const optimization [#1] INFO:Minimization -- RooAbsMinimizerFcn::setOptimizeConst: deactivating const optimization
old syntax to access the variable x RooPlot* frame = w::x.frame() ;
RooPlot *frame = x->frame();
d->plotOn(frame);
OLD syntax to omit x:: NB: The 'w::' prefix can be omitted if namespace w is imported in local namespace in the usual C++ way
using namespace w; model.plotOn(frame) ; model.plotOn(frame,Components(bkg),LineStyle(kDashed)) ;
new correct syntax
RooAbsPdf *bkg = w1->pdf("bkg");
model->plotOn(frame);
model->plotOn(frame, Components(*bkg), LineStyle(kDashed));
[#1] INFO:Plotting -- RooAbsPdf::plotOn(model) directly selected PDF components: (bkg) [#1] INFO:Plotting -- RooAbsPdf::plotOn(model) indirectly selected PDF components: ()
Draw the frame on the canvas
new TCanvas("rf509_wsinteractive", "rf509_wsinteractive", 600, 600);
gPad->SetLeftMargin(0.15);
frame->GetYaxis()->SetTitleOffset(1.4);
frame->Draw();
Draw all canvases
%jsroot on
gROOT->GetListOfCanvases()->Draw()