rf512_wsfactory_oper

Organization and simultaneous fits: operator expressions and expression-based basic pdfs in the workspace factory syntax

Author: Wouter Verkerke
This notebook tutorial was automatically generated with ROOTBOOK-izer from the macro found in the ROOT repository on Wednesday, November 30, 2022 at 11:23 AM.

In [1]:
%%cpp -d
#include "RooRealVar.h"
#include "RooDataSet.h"
#include "RooGaussian.h"
#include "RooChebychev.h"
#include "RooAddPdf.h"
#include "RooWorkspace.h"
#include "RooPlot.h"
#include "TCanvas.h"
#include "TAxis.h"
using namespace RooFit;
In [2]:
RooWorkspace *w = new RooWorkspace("w");

You can define typedefs for even shorter construction semantics

In [3]:
w->factory("$Typedef(Gaussian,Gaus)");
w->factory("$Typedef(Chebychev,Cheby)");

Operator p.d.f. examples

PDF addition is done with SUM (coef1*pdf1,pdf2)

In [4]:
w->factory("SUM::summodel( f[0,1]*Gaussian::gx(x[-10,10],m[0],1.0), Chebychev::ch(x,{0.1,0.2,-0.3}) )");

Extended PDF addition is done with SUM (yield1pdf1,yield2pdf2)

In [5]:
w->factory("SUM::extsummodel( Nsig[0,1000]*gx, Nbkg[0,1000]*ch )");

PDF multiplication is done with PROD ( pdf1, pdf2 )

In [6]:
w->factory("PROD::gxz( gx, Gaussian::gz(z[-10,10],0,1) )");

Conditional pdf multiplication is done with PROD ( pdf1|obs, pdf2 )

In [7]:
w->factory("Gaussian::gy( y[-10,10], x, 1.0 )");
w->factory("PROD::gxycond( gy|x, gx )");

Convolution (numeric/ fft) is done with NCONV/FCONV (obs,pdf1,pdf2)

In [8]:
w->factory("FCONV::lxg( x, Gaussian::g(x,mg[0],1), Landau::lc(x,0,1) )");
[#1] INFO:Caching -- Changing internal binning of variable 'x' in FFT 'lxg' from 100 to 930 to improve the precision of the numerical FFT. This can be done manually by setting an additional binning named 'cache'.

Simultaneous pdfs are constructed with SIMUL( index, state1=pdf1, state2=pdf2,...)

In [9]:
w->factory("SIMUL::smodel( c[A=0,B=1], A=Gaussian::gs(x,m,s[1]), B=Landau::ls(x,0,1) )");
[#0] WARNING:InputArguments -- The parameter 's' with range [-1e+30, 1e+30] of the RooGaussian 'gs' exceeds the safe range of (0, inf). Advise to limit its range.

Operator function examples

Function multiplication is done with prod (func1, func2,...)

In [10]:
w->factory("prod::uv(u[10],v[10])");

Function addition is done with sum(func1,func2)

In [11]:
w->factory("sum::uv2(u,v)");

Lagrangian morphing function for the example shown in rf711_lagrangianmorph

In [12]:
std::string infilename = std::string(gROOT->GetTutorialDir()) + "/roofit/input_histos_rf_lagrangianmorph.root";
w->factory(("lagrangianmorph::morph($fileName('"+infilename+"'),$observableName('pTV'),$couplings({cHq3[0,1],SM[1]}),$NewPhysics(cHq3=1,SM=0),$folders({'SM_NPsq0','cHq3_NPsq1','cHq3_NPsq2'}))").c_str());
[#0] PROGRESS:InputArguments -- initializing physics inputs from file /home/sftnight/build/workspace/root-makedoc-master/rootspi/rdoc/src/master.build/tutorials/roofit/input_histos_rf_lagrangianmorph.root with object name(s) 'pTV'
[#1] INFO:ObjectHandling -- RooWorkspace::import(w) importing dataset dh_SM_NPsq0_morph
[#1] INFO:ObjectHandling -- RooWorkspace::import(w) importing dataset dh_cHq3_NPsq1_morph
[#1] INFO:ObjectHandling -- RooWorkspace::import(w) importing dataset dh_cHq3_NPsq2_morph

Taylor expansion is done with taylorexpand(func,{var1,var2,...},val,order)

In [13]:
w->factory("taylorexpand::te(expr::poly('x^4+5*x^3+2*x^2+x+1',x),{x},0.0,2)");

Interpreted and compiled expression based p.d.f.s.

Create a RooGenericPdf interpreted pdf You can use single quotes to pass the expression string argument

In [14]:
w->factory("EXPR::G('x*x+1',x)");

Create a custom compiled pdf similar to the above interpreted pdf The code required to make this pdf is automatically embedded in the workspace

In [15]:
w->factory("CEXPR::GC('x*x+a',{x,a[1]})");
[#1] INFO:ObjectHandling -- RooWorkspace::autoImportClass(w) importing code of class RooCFAuto000Pdf from /home/sftnight/build/workspace/root-makedoc-master/rootspi/rdoc/master_TMP/notebooks/RooCFAuto000Pdf.cxx and RooCFAuto000Pdf.h
Info in <TUnixSystem::ACLiC>: creating shared library /home/sftnight/build/workspace/root-makedoc-master/rootspi/rdoc/master_TMP/notebooks/RooCFAuto000Pdf_cxx.so
In file included from input_line_9:6:
In file included from ./RooCFAuto000Pdf.cxx:11:
./RooCFAuto000Pdf.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 RooCFAuto000Pdf(*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 ./RooCFAuto000Pdf.cxx:11:
./RooCFAuto000Pdf.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;
                 ^

Compiled and interpreted functions (rather than pdfs) can be made with the lower case 'expr' and 'cexpr' types

Print workspace contents

In [16]:
w->Print();
RooWorkspace(w) w contents

variables
---------
(Nbkg,Nsig,SM,a,binWidth_pTV,c,cHq3,f,m,mg,nNP0,nNP1,nNP2,nNP3,nNP4,pTV,s,te_c0,te_c1,te_c2,te_x^0,te_x^1,te_x^2,u,v,x,y,z)

p.d.f.s
-------
RooGenericPdf::G[ actualVars=(x) formula="x[0]*x[0]+1" ] = 1
RooCFAuto000Pdf::GC[ x=x a=a ] = 1
RooChebychev::ch[ x=x coefList=(0.1,0.2,-0.3) ] = 0.8
RooAddPdf::extsummodel[ Nsig * gx + Nbkg * ch ] = 0.9/1
RooGaussian::g[ x=x mean=mg sigma=1 ] = 1
RooGaussian::gs[ x=x mean=m sigma=s ] = 1
RooGaussian::gx[ x=x mean=m sigma=1 ] = 1
RooProdPdf::gxycond[ gx * gy|x ] = 1
RooProdPdf::gxz[ gx * gz ] = 1
RooGaussian::gy[ x=y mean=x sigma=1 ] = 1
RooGaussian::gz[ x=z mean=0 sigma=1 ] = 1
RooLandau::lc[ x=x mean=0 sigma=1 ] = 0.178854
RooLandau::ls[ x=x mean=0 sigma=1 ] = 0.178854
RooFFTConvPdf::lxg[ g(x) (*) lc(x) ] = 0.171919
RooSimultaneous::smodel[ indexCat=c A=gs B=ls ] = 1
RooAddPdf::summodel[ f * gx + [%] * ch ] = 0.9/1

functions
--------
RooLagrangianMorphFunc::morph[ physics=(phys_SM_NPsq0_morph,phys_cHq3_NPsq1_morph,phys_cHq3_NPsq2_morph) operators=(cHq3,SM) observables=(pTV) binWidths=(binWidth_pTV) flags=(nNP0,nNP1,nNP2,nNP3,nNP4) binWidth_pTV * SM_NPsq0_morph + binWidth_pTV * cHq3_NPsq1_morph + binWidth_pTV * cHq3_NPsq2_morph ] = 35.9835
RooHistFunc::phys_SM_NPsq0_morph[ depList=(pTV) ] = 65.1309
RooHistFunc::phys_cHq3_NPsq1_morph[ depList=(pTV) ] = 544.966
RooHistFunc::phys_cHq3_NPsq2_morph[ depList=(pTV) ] = 1528.22
RooFormulaVar::poly[ actualVars=(x) formula="x^4+5*x^3+2*x^2+x+1" ] = 1
RooPolyFunc::te[ vars=(x) te_t0=(te_x^0,te_c0) te_t1=(te_x^1,te_c1) te_t2=(te_x^2,te_c2) ] = 1
RooProduct::uv[ u * v ] = 100
RooAddition::uv2[ u + v ] = 20

embedded datasets (in pdfs and functions)
-----------------------------------------
RooDataHist::dh_SM_NPsq0_morph(pTV)
RooDataHist::dh_cHq3_NPsq1_morph(pTV)
RooDataHist::dh_cHq3_NPsq2_morph(pTV)

embedded class code
-------------------
RooCFAuto000Pdf

embedded precalculated expensive components
-------------------------------------------
uid = 0 key=lxg_g_CONV_lc_CACHEHIST_Obs[x]_BufFrac0.1_BufStrat0 value=RooDataHist::lxg_g_CONV_lc_CACHEHIST_Obs[x]_BufFrac0.1_BufStrat0 parameters=( mg=0 )

Make workspace visible on command line

In [17]:
gDirectory->Add(w);