rf306_condpereventerrors

Multidimensional models: conditional pdf with per-event errors

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:22 AM.

In [ ]:
%%cpp -d
#include "RooRealVar.h"
#include "RooDataSet.h"
#include "RooGaussian.h"
#include "RooGaussModel.h"
#include "RooConstVar.h"
#include "RooDecay.h"
#include "RooLandau.h"
#include "RooPlot.h"
#include "TCanvas.h"
#include "TAxis.h"
#include "TH1.h"
using namespace RooFit;

B-physics pdf with per-event Gaussian resolution

Observables

In [ ]:
RooRealVar dt("dt", "dt", -10, 10);
RooRealVar dterr("dterr", "per-event error on dt", 0.01, 10);

Build a gaussian resolution model scaled by the per-event error = gauss(dt,bias,sigma*dterr)

In [ ]:
RooRealVar bias("bias", "bias", 0, -10, 10);
RooRealVar sigma("sigma", "per-event error scale factor", 1, 0.1, 10);
RooGaussModel gm("gm1", "gauss model scaled bt per-event error", dt, bias, sigma, dterr);

Construct decay(dt) (x) gauss1(dt|dterr)

In [ ]:
RooRealVar tau("tau", "tau", 1.548);
RooDecay decay_gm("decay_gm", "decay", dt, tau, gm, RooDecay::DoubleSided);

Construct fake 'external' data with per-event error

Use landau pdf to get somewhat realistic distribution with long tail

In [ ]:
RooLandau pdfDtErr("pdfDtErr", "pdfDtErr", dterr, RooConst(1), RooConst(0.25));
RooDataSet *expDataDterr = pdfDtErr.generate(dterr, 10000);

Sample data from conditional decay_gm(dt|dterr)

Specify external dataset with dterr values to use decay_dm as conditional pdf

In [ ]:
RooDataSet *data = decay_gm.generate(dt, ProtoData(*expDataDterr));

Fit conditional decay_dm(dt|dterr)

Specify dterr as conditional observable

In [ ]:
decay_gm.fitTo(*data, ConditionalObservables(dterr));

Plot conditional decay_dm(dt|dterr)

Make two-dimensional plot of conditional pdf in (dt,dterr)

In [ ]:
TH1 *hh_decay = decay_gm.createHistogram("hh_decay", dt, Binning(50), YVar(dterr, Binning(50)));
hh_decay->SetLineColor(kBlue);

Plot decay_gm(dt|dterr) at various values of dterr

In [ ]:
RooPlot *frame = dt.frame(Title("Slices of decay(dt|dterr) at various dterr"));
for (Int_t ibin = 0; ibin < 100; ibin += 20) {
   dterr.setBin(ibin);
   decay_gm.plotOn(frame, Normalization(5.));
}

Make projection of data an dt

In [ ]:
RooPlot *frame2 = dt.frame(Title("Projection of decay(dt|dterr) on dt"));
data->plotOn(frame2);

Make projection of decay(dt|dterr) on dt.

Instead of integrating out dterr, make a weighted average of curves at values dterr_i as given in the external dataset. (The true argument bins the data before projection to speed up the process)

In [ ]:
decay_gm.plotOn(frame2, ProjWData(*expDataDterr, true));

Draw all frames on canvas

In [ ]:
TCanvas *c = new TCanvas("rf306_condpereventerrors", "rf306_condperventerrors", 1200, 400);
c->Divide(3);
c->cd(1);
gPad->SetLeftMargin(0.20);
hh_decay->GetZaxis()->SetTitleOffset(2.5);
hh_decay->Draw("surf");
c->cd(2);
gPad->SetLeftMargin(0.15);
frame->GetYaxis()->SetTitleOffset(1.6);
frame->Draw();
c->cd(3);
gPad->SetLeftMargin(0.15);
frame2->GetYaxis()->SetTitleOffset(1.6);
frame2->Draw();

Draw all canvases

In [ ]:
gROOT->GetListOfCanvases()->Draw()