rf104_classfactory

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 Sunday, November 27, 2022 at 11:06 AM.

In [1]:
%%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 class skeleton code

Write skeleton pdf class with variable x,a,b To use this class,

  • Edit the file MyPdfV1.cxx and implement the evaluate() method in terms of x,a and b
  • Compile and link class with '.x MyPdfV1.cxx+'
In [2]:
RooClassFactory::makePdf("MyPdfV1", "x,A,B");

With added initial value expression

Write skeleton pdf class with variable x,a,b and given formula expression To use this class,

  • Compile and link class with '.x MyPdfV2.cxx+'
In [3]:
RooClassFactory::makePdf("MyPdfV2", "x,A,B", "", "A*fabs(x)+pow(x-B,2)");

With added analytical integral expression

Write skeleton pdf class with variable x,a,b, given formula expression and given expression for analytical integral over x To use this class,

  • Compile and link class with '.x MyPdfV3.cxx+'
In [4]:
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))");

Use instance of created class

Compile MyPdfV3 class

In [5]:
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
In file included from input_line_9:6:
In file included from ././MyPdfV3.cxx:11:
/home/sftnight/build/workspace/root-makedoc-master/rootspi/rdoc/master_TMP/notebooks/MyPdfV3.h:24:20: warning: 'clone' overrides a member function but is not marked 'override' [-Winconsistent-missing-override]
  virtual TObject* clone(const char* newname) const { return new MyPdfV3(*this,newname); }
                   ^
/home/sftnight/build/workspace/root-makedoc-master/rootspi/rdoc/src/master.build/include/RooAbsArg.h:82:20: note: overridden virtual function is here
  virtual TObject* clone(const char* newname=nullptr) const = 0 ;
                   ^
In file included from input_line_9:6:
In file included from ././MyPdfV3.cxx:11:
/home/sftnight/build/workspace/root-makedoc-master/rootspi/rdoc/master_TMP/notebooks/MyPdfV3.h:27:9: warning: 'getAnalyticalIntegral' overrides a member function but is not marked 'override' [-Winconsistent-missing-override]
  Int_t getAnalyticalIntegral(RooArgSet& allVars, RooArgSet& analVars, const char* rangeName=0) const ;
        ^
/home/sftnight/build/workspace/root-makedoc-master/rootspi/rdoc/src/master.build/include/RooAbsReal.h:166:17: note: overridden virtual function is here
  virtual Int_t getAnalyticalIntegral(RooArgSet& allVars, RooArgSet& analVars, const char* rangeName=nullptr) const ;
                ^
In file included from input_line_9:6:
In file included from ././MyPdfV3.cxx:11:
/home/sftnight/build/workspace/root-makedoc-master/rootspi/rdoc/master_TMP/notebooks/MyPdfV3.h:28:10: warning: 'analyticalIntegral' overrides a member function but is not marked 'override' [-Winconsistent-missing-override]
  double analyticalIntegral(Int_t code, const char* rangeName=0) const ;
         ^
/home/sftnight/build/workspace/root-makedoc-master/rootspi/rdoc/src/master.build/include/RooAbsReal.h:167:18: note: overridden virtual function is here
  virtual double analyticalIntegral(Int_t code, const char* rangeName=nullptr) const ;
                 ^
In file included from input_line_9:6:
In file included from ././MyPdfV3.cxx:11:
/home/sftnight/build/workspace/root-makedoc-master/rootspi/rdoc/master_TMP/notebooks/MyPdfV3.h:36:10: warning: 'evaluate' overrides a member function but is not marked 'override' [-Winconsistent-missing-override]
  double evaluate() const ;
         ^
/home/sftnight/build/workspace/root-makedoc-master/rootspi/rdoc/src/master.build/include/RooAbsReal.h:436:18: note: overridden virtual function is here
  virtual double evaluate() const = 0;
                 ^

Create instance of MyPdfV3 class

In [6]:
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

In [7]:
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

In [8]:
RooPlot *frame1 = y.frame(Title("Compiled class MyPdfV3"));
RooDataSet *data = pdf->generate(y, 1000);
pdf->fitTo(*data);
data->plotOn(frame1);
pdf->plotOn(frame1);
[#1] INFO:Minimization -- RooAbsMinimizerFcn::setOptimizeConst: activating const optimization
 **********
 **    1 **SET PRINT           1
 **********
 **********
 **    2 **SET NOGRAD
 **********
 PARAMETER DEFINITIONS:
    NO.   NAME         VALUE      STEP SIZE      LIMITS
     1 b            2.00000e+00  2.00000e+00   -1.00000e+01  1.00000e+01
 **********
 **    3 **SET ERR         0.5
 **********
 **********
 **    4 **SET PRINT           1
 **********
 **********
 **    5 **SET STR           1
 **********
 NOW USING STRATEGY  1: TRY TO BALANCE SPEED AGAINST RELIABILITY
 **********
 **    6 **MIGRAD         500           1
 **********
 FIRST CALL TO USER FUNCTION AT NEW START POINT, WITH IFLAG=4.
 START MIGRAD MINIMIZATION.  STRATEGY  1.  CONVERGENCE WHEN EDM .LT. 1.00e-03
 FCN=2568.49 FROM MIGRAD    STATUS=INITIATE        4 CALLS           5 TOTAL
                     EDM= unknown      STRATEGY= 1      NO ERROR MATRIX       
  EXT PARAMETER               CURRENT GUESS       STEP         FIRST   
  NO.   NAME      VALUE            ERROR          SIZE      DERIVATIVE 
   1  b            2.00000e+00   2.00000e+00   2.05758e-01   5.02230e+01
                               ERR DEF= 0.5
 MIGRAD MINIMIZATION HAS CONVERGED.
 MIGRAD WILL VERIFY CONVERGENCE AND ERROR MATRIX.
 COVARIANCE MATRIX CALCULATED SUCCESSFULLY
 FCN=2568.28 FROM MIGRAD    STATUS=CONVERGED      17 CALLS          18 TOTAL
                     EDM=9.63237e-09    STRATEGY= 1      ERROR MATRIX ACCURATE 
  EXT PARAMETER                                   STEP         FIRST   
  NO.   NAME      VALUE            ERROR          SIZE      DERIVATIVE 
   1  b            1.91995e+00   1.24398e-01   4.43576e-04   7.74259e-03
                               ERR DEF= 0.5
 EXTERNAL ERROR MATRIX.    NDIM=  25    NPAR=  1    ERR DEF=0.5
  1.548e-02 
 **********
 **    7 **SET ERR         0.5
 **********
 **********
 **    8 **SET PRINT           1
 **********
 **********
 **    9 **HESSE         500
 **********
 COVARIANCE MATRIX CALCULATED SUCCESSFULLY
 FCN=2568.28 FROM HESSE     STATUS=OK              5 CALLS          23 TOTAL
                     EDM=9.63457e-09    STRATEGY= 1      ERROR MATRIX ACCURATE 
  EXT PARAMETER                                INTERNAL      INTERNAL  
  NO.   NAME      VALUE            ERROR       STEP SIZE       VALUE   
   1  b            1.91995e+00   1.24398e-01   1.77431e-05   1.93194e-01
                               ERR DEF= 0.5
 EXTERNAL ERROR MATRIX.    NDIM=  25    NPAR=  1    ERR DEF=0.5
  1.548e-02 
[#1] INFO:Minimization -- RooAbsMinimizerFcn::setOptimizeConst: deactivating const optimization
input_line_64:3:1: warning: 'data' shadows a declaration with the same name in the 'std' namespace; use '::data' to reference this declaration
RooDataSet *data = pdf->generate(y, 1000);
^

Compiled version of example rf103

Declare observable x

In [9]:
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

In [10]:
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
In file included from input_line_9:6:
In file included from ./RooGenPdfPdf.cxx:11:
./RooGenPdfPdf.h:23:20: warning: 'clone' overrides a member function but is not marked 'override' [-Winconsistent-missing-override]
  virtual TObject* clone(const char* newname) const { return new RooGenPdfPdf(*this,newname); }
                   ^
/home/sftnight/build/workspace/root-makedoc-master/rootspi/rdoc/src/master.build/include/RooAbsArg.h:82:20: note: overridden virtual function is here
  virtual TObject* clone(const char* newname=nullptr) const = 0 ;
                   ^
In file included from input_line_9:6:
In file included from ./RooGenPdfPdf.cxx:11:
./RooGenPdfPdf.h:31:10: warning: 'evaluate' overrides a member function but is not marked 'override' [-Winconsistent-missing-override]
  double evaluate() const ;
         ^
/home/sftnight/build/workspace/root-makedoc-master/rootspi/rdoc/src/master.build/include/RooAbsReal.h:436:18: note: overridden virtual function is here
  virtual double evaluate() const = 0;
                 ^

Generate a toy dataset from the interpreted pdf

In [11]:
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

In [12]:
genpdf->fitTo(*data2);
[#1] INFO:NumericIntegration -- RooRealIntegral::init(GenPdf_Int[x]) using numeric integrator RooIntegrator1D to calculate Int(x)
[#1] INFO:Minimization -- RooAbsMinimizerFcn::setOptimizeConst: activating const optimization
 **********
 **   10 **SET PRINT           1
 **********
 **********
 **   11 **SET NOGRAD
 **********
 PARAMETER DEFINITIONS:
    NO.   NAME         VALUE      STEP SIZE      LIMITS
     1 alpha        5.00000e+00  9.90000e-01    1.00000e-01  1.00000e+01
 **********
 **   12 **SET ERR         0.5
 **********
 **********
 **   13 **SET PRINT           1
 **********
 **********
 **   14 **SET STR           1
 **********
 NOW USING STRATEGY  1: TRY TO BALANCE SPEED AGAINST RELIABILITY
 **********
 **   15 **MIGRAD         500           1
 **********
 FIRST CALL TO USER FUNCTION AT NEW START POINT, WITH IFLAG=4.
 START MIGRAD MINIMIZATION.  STRATEGY  1.  CONVERGENCE WHEN EDM .LT. 1.00e-03
 FCN=178267 FROM MIGRAD    STATUS=INITIATE        6 CALLS           7 TOTAL
                     EDM= unknown      STRATEGY= 1      NO ERROR MATRIX       
  EXT PARAMETER               CURRENT GUESS       STEP         FIRST   
  NO.   NAME      VALUE            ERROR          SIZE      DERIVATIVE 
   1  alpha        5.00000e+00   9.90000e-01   2.01369e-01   2.50886e+02
                               ERR DEF= 0.5
 MIGRAD MINIMIZATION HAS CONVERGED.
 MIGRAD WILL VERIFY CONVERGENCE AND ERROR MATRIX.
 COVARIANCE MATRIX CALCULATED SUCCESSFULLY
 FCN=178266 FROM MIGRAD    STATUS=CONVERGED      14 CALLS          15 TOTAL
                     EDM=3.47331e-06    STRATEGY= 1      ERROR MATRIX ACCURATE 
  EXT PARAMETER                                   STEP         FIRST   
  NO.   NAME      VALUE            ERROR          SIZE      DERIVATIVE 
   1  alpha        4.98248e+00   1.85558e-02   1.09523e-03  -4.97115e-01
                               ERR DEF= 0.5
 EXTERNAL ERROR MATRIX.    NDIM=  25    NPAR=  1    ERR DEF=0.5
  3.443e-04 
 **********
 **   16 **SET ERR         0.5
 **********
 **********
 **   17 **SET PRINT           1
 **********
 **********
 **   18 **HESSE         500
 **********
 COVARIANCE MATRIX CALCULATED SUCCESSFULLY
 FCN=178266 FROM HESSE     STATUS=OK              5 CALLS          20 TOTAL
                     EDM=3.46509e-06    STRATEGY= 1      ERROR MATRIX ACCURATE 
  EXT PARAMETER                                INTERNAL      INTERNAL  
  NO.   NAME      VALUE            ERROR       STEP SIZE       VALUE   
   1  alpha        4.98248e+00   1.85557e-02   2.19047e-04  -1.36416e-02
                               ERR DEF= 0.5
 EXTERNAL ERROR MATRIX.    NDIM=  25    NPAR=  1    ERR DEF=0.5
  3.443e-04 
[#1] INFO:Minimization -- RooAbsMinimizerFcn::setOptimizeConst: deactivating const optimization

Make a plot of the data and the pdf overlaid

In [13]:
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

In [14]:
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

In [15]:
%jsroot on
gROOT->GetListOfCanvases()->Draw()