Basic functionality: The class factory for functions and pdfs
NOTE: This demo uses code that is generated by the macro, therefore it cannot be compiled in one step by ACliC. To run this macro compiled with ACliC do
root>.x rf104_classfactory.C // run interpreted to generate code
root>.L MyPdfV3.cxx+ // Compile and load created class
root>.x rf104_classfactory.C+ // run compiled code
Author: Wouter Verkerke
This notebook tutorial was automatically generated with ROOTBOOK-izer from the macro found in the ROOT repository on Wednesday, April 17, 2024 at 11:16 AM.
%%cpp -d
#include "RooRealVar.h"
#include "RooDataSet.h"
#include "RooGaussian.h"
#include "TCanvas.h"
#include "TAxis.h"
#include "RooPlot.h"
#include "RooClassFactory.h"
#include "TROOT.h"
using namespace RooFit;
Write skeleton pdf class with variable x,a,b To use this class,
RooClassFactory::makePdf("MyPdfV1", "x,A,B");
Write skeleton pdf class with variable x,a,b and given formula expression To use this class,
RooClassFactory::makePdf("MyPdfV2", "x,A,B", "", "A*fabs(x)+pow(x-B,2)");
Write skeleton pdf class with variable x,a,b, given formula expression and given expression for analytical integral over x To use this class,
RooClassFactory::makePdf("MyPdfV3", "x,A,B", "", "A*fabs(x)+pow(x-B,2)", true, false,
"x:(A/2)*(pow(x.max(rangeName),2)+pow(x.min(rangeName),2))+(1./"
"3)*(pow(x.max(rangeName)-B,3)-pow(x.min(rangeName)-B,3))");
Compile MyPdfV3 class
gROOT->ProcessLineSync(".x MyPdfV3.cxx+");
(MyPdfV3) An instance of MyPdfV3.
Info in <TUnixSystem::ACLiC>: creating shared library /home/sftnight/build/workspace/root-makedoc-master/rootspi/rdoc/master_TMP/notebooks/./MyPdfV3_cxx.so c++: error: /home/sftnight/build/workspace/root-makedoc-master/rootspi/rdoc/master_TMP/notebooks/MyPdfV3_cxx_ACLiC_dict.cxx: No such file or directory c++: fatal error: no input files compilation terminated.
Create instance of MyPdfV3 class
RooRealVar a("a", "a", 1);
RooRealVar b("b", "b", 2, -10, 10);
RooRealVar y("y", "y", -10, 10);
We need to hide the type to run in a ROOT macro
RooWorkspace w("w");
w.factory("MyPdfV3::pdf(y[-10,10], a[1], b[2,-10,10])");
auto pdf = w.pdf("pdf");
Generate toy data from pdf and plot data and pdf on frame
RooPlot *frame1 = y.frame(Title("Compiled class MyPdfV3"));
std::unique_ptr<RooDataSet> data{pdf->generate(y, 1000)};
pdf->fitTo(*data, PrintLevel(-1));
data->plotOn(frame1);
pdf->plotOn(frame1);
[#1] INFO:Fitting -- RooAbsPdf::fitTo(pdf_over_pdf_Int[y]) 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_pdf_over_pdf_Int[y]_pdfData) Summation contains a RooNLLVar, using its error level [#1] INFO:Minimization -- RooAbsMinimizerFcn::setOptimizeConst: activating const optimization [#1] INFO:Minimization -- RooAbsMinimizerFcn::setOptimizeConst: deactivating const optimization
input_line_63:3:1: warning: 'data' shadows a declaration with the same name in the 'std' namespace; use '::data' to reference this declaration std::unique_ptr<RooDataSet> data{pdf->generate(y, 1000)}; ^
Declare observable x
RooRealVar x("x", "x", -20, 20);
The RooClassFactory::makePdfInstance() function performs code writing, compiling, linking and object instantiation in one go and can serve as a straight replacement of RooGenericPdf
RooRealVar alpha("alpha", "alpha", 5, 0.1, 10);
RooAbsPdf *genpdf =
RooClassFactory::makePdfInstance("GenPdf", "(1+0.1*fabs(x)+sin(sqrt(fabs(x*alpha+0.1))))", RooArgSet(x, alpha));
Info in <TUnixSystem::ACLiC>: creating shared library /home/sftnight/build/workspace/root-makedoc-master/rootspi/rdoc/master_TMP/notebooks/RooGenPdfPdf_cxx.so
Generate a toy dataset from the interpreted pdf
std::unique_ptr<RooDataSet> data2{genpdf->generate(x, 50000)};
[#1] INFO:NumericIntegration -- RooRealIntegral::init(GenPdf_Int[x]) using numeric integrator RooIntegrator1D to calculate Int(x) [#1] INFO:NumericIntegration -- RooRealIntegral::init(GenPdf_Int[x]) using numeric integrator RooIntegrator1D to calculate Int(x)
Fit the interpreted pdf to the generated data
genpdf->fitTo(*data2, PrintLevel(-1));
[#1] INFO:Fitting -- RooAbsPdf::fitTo(GenPdf_over_GenPdf_Int[x]) fixing normalization set for coefficient determination to observables in data [#1] INFO:Fitting -- RooAddition::defaultErrorLevel(nll_GenPdf_over_GenPdf_Int[x]_GenPdfData) Summation contains a RooNLLVar, using its error level [#1] INFO:Minimization -- RooAbsMinimizerFcn::setOptimizeConst: activating const optimization [#1] INFO:NumericIntegration -- RooRealIntegral::init(GenPdf_Int[x]) using numeric integrator RooIntegrator1D to calculate Int(x) [#1] INFO:Minimization -- RooAbsMinimizerFcn::setOptimizeConst: deactivating const optimization
Make a plot of the data and the pdf overlaid
RooPlot *frame2 = x.frame(Title("Compiled version of pdf of rf103"));
data2->plotOn(frame2);
genpdf->plotOn(frame2);
[#1] INFO:NumericIntegration -- RooRealIntegral::init(GenPdf_Int[x]) using numeric integrator RooIntegrator1D to calculate Int(x)
Draw all frames on a canvas
TCanvas *c = new TCanvas("rf104_classfactory", "rf104_classfactory", 800, 400);
c->Divide(2);
c->cd(1);
gPad->SetLeftMargin(0.15);
frame1->GetYaxis()->SetTitleOffset(1.4);
frame1->Draw();
c->cd(2);
gPad->SetLeftMargin(0.15);
frame2->GetYaxis()->SetTitleOffset(1.4);
frame2->Draw();
Draw all canvases
%jsroot on
gROOT->GetListOfCanvases()->Draw()