Logo ROOT  
Reference Guide
rf210_angularconv.C File Reference

Detailed Description

View in nbviewer Open in SWAN

Addition and convolution: convolution in cyclical angular observables theta

and construction of p.d.f in terms of transformed angular coordinates, e.g. cos(theta), where the convolution is performed in theta rather than cos(theta)

pdf(theta) = T(theta) (x) gauss(theta)
pdf(cosTheta) = T(acos(cosTheta)) (x) gauss(acos(cosTheta))

This tutorial requires FFT3 to be enabled.

RooFit v3.60 -- Developed by Wouter Verkerke and David Kirkby
Copyright (C) 2000-2013 NIKHEF, University of California & Stanford University
All rights reserved, please read http://roofit.sourceforge.net/license.txt
[#1] INFO:Caching -- Changing internal binning of variable 'psi' in FFT 'Mf' from 100 to 930 to improve the precision of the numerical FFT. This can be done manually by setting an additional binning named 'cache'.
[#1] INFO:Eval -- RooRealVar::setRange(psi) new range named 'refrange_fft_Mf' created with bounds [0,3.14159]
[#1] INFO:NumericIntegration -- RooRealIntegral::init(Tpsi_Int[psi]) using numeric integrator RooIntegrator1D to calculate Int(psi)
[#1] INFO:Caching -- RooAbsCachedPdf::getCache(Mf) creating new cache 0x56107a842590 with pdf Tpsi_CONV_Rpsi_CACHE_Obs[psi] for nset (psi) with code 0
[#1] INFO:Caching -- RooAbsCachedPdf::getCache(Mf) creating new cache 0x56107c696610 with pdf Tpsi_CONV_Rpsi_CACHE_Obs[psi] for nset (psi) with code 0 from preexisting content.
[#1] INFO:Caching -- RooAbsCachedPdf::getCache(Mf) creating new cache 0x56107a842590 with pdf Tpsi_CONV_Rpsi_CACHE_Obs[psi] for nset (psi) with code 0 from preexisting content.
[#1] INFO:Minization -- RooMinimizer::optimizeConst: activating const optimization
[#1] INFO:NumericIntegration -- RooRealIntegral::init(Tpsi_Int[psi]) using numeric integrator RooIntegrator1D to calculate Int(psi)
[#1] INFO:Minization -- The following expressions have been identified as constant and will be precalculated and cached: (Tpsi)
**********
** 1 **SET PRINT 1
**********
**********
** 2 **SET NOGRAD
**********
PARAMETER DEFINITIONS:
NO. NAME VALUE STEP SIZE LIMITS
1 gbias 2.00000e-01 1.00000e-01 0.00000e+00 1.00000e+00
2 greso 3.00000e-01 9.00000e-02 1.00000e-01 1.00000e+00
**********
** 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 1000 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
[#1] INFO:NumericIntegration -- RooRealIntegral::init(Tpsi_Int[psi]) using numeric integrator RooIntegrator1D to calculate Int(psi)
FCN=9480.16 FROM MIGRAD STATUS=INITIATE 8 CALLS 9 TOTAL
EDM= unknown STRATEGY= 1 NO ERROR MATRIX
EXT PARAMETER CURRENT GUESS STEP FIRST
NO. NAME VALUE ERROR SIZE DERIVATIVE
1 gbias 2.00000e-01 1.00000e-01 2.57889e-01 5.40962e+01
2 greso 3.00000e-01 9.00000e-02 2.46497e-01 5.76705e+00
ERR DEF= 0.5
MIGRAD MINIMIZATION HAS CONVERGED.
MIGRAD WILL VERIFY CONVERGENCE AND ERROR MATRIX.
COVARIANCE MATRIX CALCULATED SUCCESSFULLY
FCN=9479.64 FROM MIGRAD STATUS=CONVERGED 32 CALLS 33 TOTAL
EDM=1.07647e-07 STRATEGY= 1 ERROR MATRIX ACCURATE
EXT PARAMETER STEP FIRST
NO. NAME VALUE ERROR SIZE DERIVATIVE
1 gbias 1.92485e-01 7.44189e-03 1.26877e-03 1.36618e-02
2 greso 2.98582e-01 9.19693e-03 1.65656e-03 -8.25561e-03
ERR DEF= 0.5
EXTERNAL ERROR MATRIX. NDIM= 25 NPAR= 2 ERR DEF=0.5
5.539e-05 1.690e-07
1.690e-07 8.460e-05
PARAMETER CORRELATION COEFFICIENTS
NO. GLOBAL 1 2
1 0.00247 1.000 0.002
2 0.00247 0.002 1.000
**********
** 7 **SET ERR 0.5
**********
**********
** 8 **SET PRINT 1
**********
**********
** 9 **HESSE 1000
**********
COVARIANCE MATRIX CALCULATED SUCCESSFULLY
FCN=9479.64 FROM HESSE STATUS=OK 10 CALLS 43 TOTAL
EDM=1.07707e-07 STRATEGY= 1 ERROR MATRIX ACCURATE
EXT PARAMETER INTERNAL INTERNAL
NO. NAME VALUE ERROR STEP SIZE VALUE
1 gbias 1.92485e-01 7.44188e-03 5.07510e-05 -6.62425e-01
2 greso 2.98582e-01 9.19691e-03 6.62625e-05 -5.92825e-01
ERR DEF= 0.5
EXTERNAL ERROR MATRIX. NDIM= 25 NPAR= 2 ERR DEF=0.5
5.539e-05 1.248e-07
1.248e-07 8.460e-05
PARAMETER CORRELATION COEFFICIENTS
NO. GLOBAL 1 2
1 0.00182 1.000 0.002
2 0.00182 0.002 1.000
[#1] INFO:Minization -- RooMinimizer::optimizeConst: deactivating const optimization
[#1] INFO:NumericIntegration -- RooRealIntegral::init(Tpsi_Int[psi]) using numeric integrator RooIntegrator1D to calculate Int(psi)
[#1] INFO:Caching -- RooAbsCachedPdf::getCache(Mf) creating new cache 0x56107a842590 with pdf Tpsi_CONV_Rpsi_CACHE_Obs[psi] for nset (psi) with code 0
[#1] INFO:NumericIntegration -- RooRealIntegral::init(Tpsi_Int[psi]) using numeric integrator RooIntegrator1D to calculate Int(psi)
[#1] INFO:Caching -- RooAbsCachedPdf::getCache(Mf) creating new cache 0x56107c85d770 with pdf Tpsi_CONV_Rpsi_CACHE_Obs[cpsi] for nset (cpsi) with code 0 from preexisting content.
[#1] INFO:NumericIntegration -- RooRealIntegral::init(Tpsi_CONV_Rpsi_CACHE_Obs[cpsi]_Int[cpsi]) using numeric integrator RooIntegrator1D to calculate Int(cpsi)
[#1] INFO:Caching -- RooAbsCachedPdf::getCache(Mf) creating new cache 0x56107c85d680 with pdf Tpsi_CONV_Rpsi_CACHE_Obs[cpsi] for nset (cpsi) with code 0 from preexisting content.
[#1] INFO:NumericIntegration -- RooRealIntegral::init(Tpsi_CONV_Rpsi_CACHE_Obs[cpsi]_Int[cpsi]) using numeric integrator RooIntegrator1D to calculate Int(cpsi)
[#1] INFO:Caching -- RooAbsCachedPdf::getCache(Mf) creating new cache 0x56107c85d680 with pdf Tpsi_CONV_Rpsi_CACHE_Obs[cpsi] for nset (cpsi) with code 0 from preexisting content.
[#1] INFO:NumericIntegration -- RooRealIntegral::init(Tpsi_CONV_Rpsi_CACHE_Obs[cpsi]_Int[cpsi]) using numeric integrator RooIntegrator1D to calculate Int(cpsi)
[#1] INFO:Minization -- RooMinimizer::optimizeConst: activating const optimization
[#1] INFO:Minization -- The following expressions have been identified as constant and will be precalculated and cached: (psif)
**********
** 10 **SET PRINT 1
**********
**********
** 11 **SET NOGRAD
**********
PARAMETER DEFINITIONS:
NO. NAME VALUE STEP SIZE LIMITS
1 gbias 1.92485e-01 7.44188e-03 0.00000e+00 1.00000e+00
2 greso 2.98582e-01 9.19691e-03 1.00000e-01 1.00000e+00
**********
** 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 1000 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
[#1] INFO:NumericIntegration -- RooRealIntegral::init(Tpsi_Int[psi]) using numeric integrator RooIntegrator1D to calculate Int(psi)
FCN=5178.91 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 gbias 1.92485e-01 7.44188e-03 1.88791e-02 2.79436e+01
2 greso 2.98582e-01 9.19691e-03 2.46483e-02 -8.41273e+00
ERR DEF= 0.5
MIGRAD MINIMIZATION HAS CONVERGED.
MIGRAD WILL VERIFY CONVERGENCE AND ERROR MATRIX.
COVARIANCE MATRIX CALCULATED SUCCESSFULLY
FCN=5178.63 FROM MIGRAD STATUS=CONVERGED 28 CALLS 29 TOTAL
EDM=3.37796e-07 STRATEGY= 1 ERROR MATRIX ACCURATE
EXT PARAMETER STEP FIRST
NO. NAME VALUE ERROR SIZE DERIVATIVE
1 gbias 1.85814e-01 9.36461e-03 1.16493e-03 2.44691e-02
2 greso 3.02469e-01 1.02714e-02 1.32671e-03 1.52048e-03
ERR DEF= 0.5
EXTERNAL ERROR MATRIX. NDIM= 25 NPAR= 2 ERR DEF=0.5
8.771e-05 -2.174e-05
-2.174e-05 1.055e-04
PARAMETER CORRELATION COEFFICIENTS
NO. GLOBAL 1 2
1 0.22594 1.000 -0.226
2 0.22594 -0.226 1.000
**********
** 16 **SET ERR 0.5
**********
**********
** 17 **SET PRINT 1
**********
**********
** 18 **HESSE 1000
**********
COVARIANCE MATRIX CALCULATED SUCCESSFULLY
FCN=5178.63 FROM HESSE STATUS=OK 10 CALLS 39 TOTAL
EDM=3.36905e-07 STRATEGY= 1 ERROR MATRIX ACCURATE
EXT PARAMETER INTERNAL INTERNAL
NO. NAME VALUE ERROR STEP SIZE VALUE
1 gbias 1.85814e-01 9.36561e-03 2.32986e-04 -6.79459e-01
2 greso 3.02469e-01 1.02725e-02 5.30684e-05 -5.82446e-01
ERR DEF= 0.5
EXTERNAL ERROR MATRIX. NDIM= 25 NPAR= 2 ERR DEF=0.5
8.773e-05 -2.179e-05
-2.179e-05 1.056e-04
PARAMETER CORRELATION COEFFICIENTS
NO. GLOBAL 1 2
1 0.22639 1.000 -0.226
2 0.22639 -0.226 1.000
[#1] INFO:Minization -- RooMinimizer::optimizeConst: deactivating const optimization
[#1] INFO:NumericIntegration -- RooRealIntegral::init(Tpsi_Int[psi]) using numeric integrator RooIntegrator1D to calculate Int(psi)
[#1] INFO:Caching -- RooAbsCachedPdf::getCache(Mf) creating new cache 0x56107c85d680 with pdf Tpsi_CONV_Rpsi_CACHE_Obs[cpsi] for nset (cpsi) with code 0
[#1] INFO:NumericIntegration -- RooRealIntegral::init(Tpsi_CONV_Rpsi_CACHE_Obs[cpsi]_Int[cpsi]) using numeric integrator RooIntegrator1D to calculate Int(cpsi)
[#1] INFO:NumericIntegration -- RooRealIntegral::init(T_Int[cpsi]) using numeric integrator RooIntegrator1D to calculate Int(cpsi)
#include "RooRealVar.h"
#include "RooDataSet.h"
#include "RooGaussian.h"
#include "RooGenericPdf.h"
#include "RooFormulaVar.h"
#include "RooFFTConvPdf.h"
#include "RooPlot.h"
#include "TCanvas.h"
#include "TAxis.h"
#include "TH1.h"
using namespace RooFit;
void rf210_angularconv()
{
// S e t u p c o m p o n e n t p d f s
// ---------------------------------------
// Define angle psi
RooRealVar psi("psi", "psi", 0, 3.14159268);
// Define physics p.d.f T(psi)
RooGenericPdf Tpsi("Tpsi", "1+sin(2*@0)", psi);
// Define resolution R(psi)
RooRealVar gbias("gbias", "gbias", 0.2, 0., 1);
RooRealVar greso("greso", "greso", 0.3, 0.1, 1.0);
RooGaussian Rpsi("Rpsi", "Rpsi", psi, gbias, greso);
// Define cos(psi) and function psif that calculates psi from cos(psi)
RooRealVar cpsi("cpsi", "cos(psi)", -1, 1);
RooFormulaVar psif("psif", "acos(cpsi)", cpsi);
// Define physics p.d.f. also as function of cos(psi): T(psif(cpsi)) = T(cpsi) ;
RooGenericPdf Tcpsi("T", "1+sin(2*@0)", psif);
// C o n s t r u c t c o n v o l u t i o n p d f i n p s i
// --------------------------------------------------------------
// Define convoluted p.d.f. as function of psi: M=[T(x)R](psi) = M(psi)
RooFFTConvPdf Mpsi("Mf", "Mf", psi, Tpsi, Rpsi);
// Set the buffer fraction to zero to obtain a true cyclical convolution
Mpsi.setBufferFraction(0);
// S a m p l e , f i t a n d p l o t c o n v o l u t e d p d f ( p s i )
// --------------------------------------------------------------------------------
// Generate some events in observable psi
RooDataSet *data_psi = Mpsi.generate(psi, 10000);
// Fit convoluted model as function of angle psi
Mpsi.fitTo(*data_psi);
// Plot cos(psi) frame with Mf(cpsi)
RooPlot *frame1 = psi.frame(Title("Cyclical convolution in angle psi"));
data_psi->plotOn(frame1);
Mpsi.plotOn(frame1);
// Overlay comparison to unsmeared physics p.d.f T(psi)
Tpsi.plotOn(frame1, LineColor(kRed));
// C o n s t r u c t c o n v o l u t i o n p d f i n c o s ( p s i )
// --------------------------------------------------------------------------
// Define convoluted p.d.f. as function of cos(psi): M=[T(x)R](psif(cpsi)) = M(cpsi)
//
// Need to give both observable psi here (for definition of convolution)
// and function psif here (for definition of observables, ultimately in cpsi)
RooFFTConvPdf Mcpsi("Mf", "Mf", psif, psi, Tpsi, Rpsi);
// Set the buffer fraction to zero to obtain a true cyclical convolution
Mcpsi.setBufferFraction(0);
// S a m p l e , f i t a n d p l o t c o n v o l u t e d p d f ( c o s p s i )
// --------------------------------------------------------------------------------
// Generate some events
RooDataSet *data_cpsi = Mcpsi.generate(cpsi, 10000);
// set psi constant to exclude to be a parameter of the fit
psi.setConstant(true);
// Fit convoluted model as function of cos(psi)
Mcpsi.fitTo(*data_cpsi);
// Plot cos(psi) frame with Mf(cpsi)
RooPlot *frame2 = cpsi.frame(Title("Same convolution in psi, expressed in cos(psi)"));
data_cpsi->plotOn(frame2);
Mcpsi.plotOn(frame2);
// Overlay comparison to unsmeared physics p.d.f Tf(cpsi)
Tcpsi.plotOn(frame2, LineColor(kRed));
// Draw frame on canvas
TCanvas *c = new TCanvas("rf210_angularconv", "rf210_angularconv", 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();
}
Date
04/2009
Author
Wouter Verkerke

Definition in file rf210_angularconv.C.

RooFormulaVar.h
c
#define c(i)
Definition: RSha256.hxx:119
RooPlot::Draw
virtual void Draw(Option_t *options=0)
Draw this plot and all of the elements it contains.
Definition: RooPlot.cxx:691
RooGaussian.h
x
Double_t x[n]
Definition: legend1.C:17
RooGaussian
Definition: RooGaussian.h:25
TCanvas.h
RooGenericPdf.h
RooDataSet.h
RooPlot::frame
static RooPlot * frame(const RooAbsRealLValue &var, Double_t xmin, Double_t xmax, Int_t nBins)
Create a new frame for a given variable in x.
Definition: RooPlot.cxx:249
RooFormulaVar
Definition: RooFormulaVar.h:30
TGeant4Unit::gauss
static constexpr double gauss
Definition: TGeant4SystemOfUnits.h:269
RooFit
Definition: RooCFunction1Binding.h:29
RooAbsData::plotOn
virtual RooPlot * plotOn(RooPlot *frame, const RooCmdArg &arg1=RooCmdArg::none(), const RooCmdArg &arg2=RooCmdArg::none(), const RooCmdArg &arg3=RooCmdArg::none(), const RooCmdArg &arg4=RooCmdArg::none(), const RooCmdArg &arg5=RooCmdArg::none(), const RooCmdArg &arg6=RooCmdArg::none(), const RooCmdArg &arg7=RooCmdArg::none(), const RooCmdArg &arg8=RooCmdArg::none()) const
Definition: RooAbsData.cxx:547
RooPlot.h
RooPlot::GetYaxis
TAxis * GetYaxis() const
Definition: RooPlot.cxx:1256
RooPlot
Definition: RooPlot.h:44
RooRealVar.h
kRed
@ kRed
Definition: Rtypes.h:66
acos
double acos(double)
RooFit::LineColor
RooCmdArg LineColor(Color_t color)
Definition: RooGlobalFunc.cxx:56
RooFFTConvPdf
PDF for the numerical (FFT) convolution of two PDFs.
Definition: RooFFTConvPdf.h:25
TCanvas
Definition: TCanvas.h:23
TAxis.h
RooGenericPdf
Definition: RooGenericPdf.h:25
ROOT::Math::Chebyshev::T
double T(double x)
Definition: ChebyshevPol.h:52
RooFFTConvPdf.h
gPad
#define gPad
Definition: TVirtualPad.h:287
RooDataSet
Definition: RooDataSet.h:33
TAttAxis::SetTitleOffset
virtual void SetTitleOffset(Float_t offset=1)
Set distance between the axis and the axis title.
Definition: TAttAxis.cxx:293
RooRealVar
Definition: RooRealVar.h:35
TH1.h
RooFit::Title
RooCmdArg Title(const char *name)
Definition: RooGlobalFunc.cxx:173