Skip to content

Commit

Permalink
[RF] Use TRandom3 explicitly in RooFit
Browse files Browse the repository at this point in the history
This demonstrates the good practice of not relying of global variables.
  • Loading branch information
guitargeek committed Nov 25, 2024
1 parent 491b86b commit ab46cb8
Show file tree
Hide file tree
Showing 8 changed files with 62 additions and 46 deletions.
23 changes: 13 additions & 10 deletions tutorials/roofit/rf102_dataimport.C
Original file line number Diff line number Diff line change
Expand Up @@ -19,22 +19,25 @@
#include "TTree.h"
#include "TH1D.h"
#include "TRandom.h"
using namespace RooFit;

TH1 *makeTH1();
TTree *makeTTree();
TH1 *makeTH1(TRandom &trnd);
TTree *makeTTree(TRandom &trnd);

void rf102_dataimport()
{
using namespace RooFit;

// ---------------------------------------------------
// I m p o r t i n g R O O T h i s t o g r a m s
// ===================================================

TRandom3 trnd{};

// I m p o r t T H 1 i n t o a R o o D a t a H i s t
// ---------------------------------------------------------

// Create a ROOT TH1 histogram
TH1 *hh = makeTH1();
TH1 *hh = makeTH1(trnd);

// Declare observable x
RooRealVar x("x", "x", -10, 10);
Expand Down Expand Up @@ -80,7 +83,7 @@ void rf102_dataimport()
// I m p o r t T T r e e i n t o a R o o D a t a S e t
// -----------------------------------------------------------

TTree *tree = makeTTree();
TTree *tree = makeTTree(trnd);

// Define 2nd observable y
RooRealVar y("y", "y", -10, 10);
Expand Down Expand Up @@ -166,26 +169,26 @@ void rf102_dataimport()
}

// Create ROOT TH1 filled with a Gaussian distribution
TH1 *makeTH1()
TH1 *makeTH1(TRandom &trnd)
{
TH1D *hh = new TH1D("hh", "hh", 25, -10, 10);
for (int i = 0; i < 100; i++) {
hh->Fill(gRandom->Gaus(0, 3));
hh->Fill(trnd.Gaus(0, 3));
}
return hh;
}

// Create ROOT TTree filled with a Gaussian distribution in x and a uniform distribution in y
TTree *makeTTree()
TTree *makeTTree(TRandom &trnd)
{
TTree *tree = new TTree("tree", "tree");
double *px = new double;
double *py = new double;
tree->Branch("x", px, "x/D");
tree->Branch("y", py, "y/D");
for (int i = 0; i < 100; i++) {
*px = gRandom->Gaus(0, 3);
*py = gRandom->Uniform() * 30 - 15;
*px = trnd.Gaus(0, 3);
*py = trnd.Uniform() * 30 - 15;
tree->Fill();
}
return tree;
Expand Down
15 changes: 8 additions & 7 deletions tutorials/roofit/rf102_dataimport.py
Original file line number Diff line number Diff line change
Expand Up @@ -15,17 +15,17 @@
from array import array


def makeTH1():
def makeTH1(trnd):

# Create ROOT ROOT.TH1 filled with a Gaussian distribution

hh = ROOT.TH1D("hh", "hh", 25, -10, 10)
for i in range(100):
hh.Fill(ROOT.gRandom.Gaus(0, 3))
hh.Fill(trnd.Gaus(0, 3))
return hh


def makeTTree():
def makeTTree(trnd):
# Create ROOT ROOT.TTree filled with a Gaussian distribution in x and a
# uniform distribution in y

Expand All @@ -35,19 +35,20 @@ def makeTTree():
tree.Branch("x", px, "x/D")
tree.Branch("y", py, "y/D")
for i in range(100):
px[0] = ROOT.gRandom.Gaus(0, 3)
py[0] = ROOT.gRandom.Uniform() * 30 - 15
px[0] = trnd.Gaus(0, 3)
py[0] = trnd.Uniform() * 30 - 15
tree.Fill()
return tree

trnd = ROOT.TRandom3()

############################
# Importing ROOT histograms
############################
# Import ROOT TH1 into a RooDataHist
# ---------------------------------------------------------
# Create a ROOT TH1 histogram
hh = makeTH1()
hh = makeTH1(trnd)

# Declare observable x
x = ROOT.RooRealVar("x", "x", -10, 10)
Expand Down Expand Up @@ -90,7 +91,7 @@ def makeTTree():
# -----------------------------------------------------------
# Import ROOT TTree into a RooDataSet

tree = makeTTree()
tree = makeTTree(trnd)

# Define 2nd observable y
y = ROOT.RooRealVar("y", "y", -10, 10)
Expand Down
6 changes: 4 additions & 2 deletions tutorials/roofit/rf303_conditional.C
Original file line number Diff line number Diff line change
Expand Up @@ -91,15 +91,17 @@ void rf303_conditional()

RooDataSet *makeFakeDataXY()
{
TRandom3 trnd{};

RooRealVar x("x", "x", -10, 10);
RooRealVar y("y", "y", -10, 10);
RooArgSet coord(x, y);

RooDataSet *d = new RooDataSet("d", "d", RooArgSet(x, y));

for (int i = 0; i < 10000; i++) {
double tmpy = gRandom->Gaus(0, 10);
double tmpx = gRandom->Gaus(0.5 * tmpy, 1);
double tmpy = trnd.Gaus(0, 10);
double tmpx = trnd.Gaus(0.5 * tmpy, 1);
if (fabs(tmpy) < 10 && fabs(tmpx) < 10) {
x.setVal(tmpx);
y.setVal(tmpy);
Expand Down
7 changes: 5 additions & 2 deletions tutorials/roofit/rf303_conditional.py
Original file line number Diff line number Diff line change
Expand Up @@ -17,15 +17,18 @@


def makeFakeDataXY():

trnd = ROOT.TRandom3()

x = ROOT.RooRealVar("x", "x", -10, 10)
y = ROOT.RooRealVar("y", "y", -10, 10)
coord = {x, y}

d = ROOT.RooDataSet("d", "d", coord)

for i in range(10000):
tmpy = ROOT.gRandom.Gaus(0, 10)
tmpx = ROOT.gRandom.Gaus(0.5 * tmpy, 1)
tmpy = trnd.Gaus(0, 10)
tmpx = trnd.Gaus(0.5 * tmpy, 1)
if (abs(tmpy) < 10) and (abs(tmpx) < 10):
x.setVal(tmpx)
y.setVal(tmpy)
Expand Down
26 changes: 14 additions & 12 deletions tutorials/roofit/rf401_importttreethx.C
Original file line number Diff line number Diff line change
Expand Up @@ -26,18 +26,20 @@

using namespace RooFit;

TH1 *makeTH1(const char *name, double mean, double sigma);
TTree *makeTTree();
TH1 *makeTH1(TRandom &trnd, const char *name, double mean, double sigma);
TTree *makeTTree(TRandom &trnd);

void rf401_importttreethx()
{
TRandom3 trnd{};

// I m p o r t m u l t i p l e T H 1 i n t o a R o o D a t a H i s t
// --------------------------------------------------------------------------

// Create thee ROOT TH1 histograms
TH1 *hh_1 = makeTH1("hh1", 0, 3);
TH1 *hh_2 = makeTH1("hh2", -3, 1);
TH1 *hh_3 = makeTH1("hh3", +3, 4);
TH1 *hh_1 = makeTH1(trnd, "hh1", 0, 3);
TH1 *hh_2 = makeTH1(trnd, "hh2", -3, 1);
TH1 *hh_3 = makeTH1(trnd, "hh3", +3, 4);

// Declare observable x
RooRealVar x("x", "x", -10, 10);
Expand All @@ -61,7 +63,7 @@ void rf401_importttreethx()
// I m p o r t i n g a T T r e e i n t o a R o o D a t a S e t w i t h c u t s
// -----------------------------------------------------------------------------------------

TTree *tree = makeTTree();
TTree *tree = makeTTree(trnd);

// Define observables y,z
RooRealVar y("y", "y", -10, 10);
Expand Down Expand Up @@ -108,18 +110,18 @@ void rf401_importttreethx()
dsABC.Print();
}

TH1 *makeTH1(const char *name, double mean, double sigma)
TH1 *makeTH1(TRandom &trnd, const char *name, double mean, double sigma)
{
// Create ROOT TH1 filled with a Gaussian distribution

TH1D *hh = new TH1D(name, name, 100, -10, 10);
for (int i = 0; i < 1000; i++) {
hh->Fill(gRandom->Gaus(mean, sigma));
hh->Fill(trnd.Gaus(mean, sigma));
}
return hh;
}

TTree *makeTTree()
TTree *makeTTree(TRandom &trnd)
{
// Create ROOT TTree filled with a Gaussian distribution in x and a uniform distribution in y

Expand All @@ -133,9 +135,9 @@ TTree *makeTTree()
tree->Branch("z", pz, "z/D");
tree->Branch("i", pi, "i/I");
for (int i = 0; i < 100; i++) {
*px = gRandom->Gaus(0, 3);
*py = gRandom->Uniform() * 30 - 15;
*pz = gRandom->Gaus(0, 5);
*px = trnd.Gaus(0, 3);
*py = trnd.Uniform() * 30 - 15;
*pz = trnd.Gaus(0, 5);
*pi = i % 3;
tree->Fill();
}
Expand Down
21 changes: 11 additions & 10 deletions tutorials/roofit/rf401_importttreethx.py
Original file line number Diff line number Diff line change
Expand Up @@ -16,17 +16,17 @@
from array import array


def makeTH1(name, mean, sigma):
def makeTH1(trnd, name, mean, sigma):
"""Create ROOT TH1 filled with a Gaussian distribution."""

hh = ROOT.TH1D(name, name, 100, -10, 10)
for i in range(1000):
hh.Fill(ROOT.gRandom.Gaus(mean, sigma))
hh.Fill(trnd.Gaus(mean, sigma))

return hh


def makeTTree():
def makeTTree(trnd):
"""Create ROOT ROOT.TTree filled with a Gaussian distribution in x and a uniform distribution in y."""

tree = ROOT.TTree("tree", "tree")
Expand All @@ -39,22 +39,23 @@ def makeTTree():
tree.Branch("z", pz, "z/D")
tree.Branch("i", pi, "i/I")
for i in range(100):
px[0] = ROOT.gRandom.Gaus(0, 3)
py[0] = ROOT.gRandom.Uniform() * 30 - 15
pz[0] = ROOT.gRandom.Gaus(0, 5)
px[0] = trnd.Gaus(0, 3)
py[0] = trnd.Uniform() * 30 - 15
pz[0] = trnd.Gaus(0, 5)
pi[0] = i % 3
tree.Fill()

return tree

trnd = ROOT.TRandom3()

# Import multiple TH1 into a RooDataHist
# ----------------------------------------------------------

# Create thee ROOT ROOT.TH1 histograms
hh_1 = makeTH1("hh1", 0, 3)
hh_2 = makeTH1("hh2", -3, 1)
hh_3 = makeTH1("hh3", +3, 4)
hh_1 = makeTH1(trnd, "hh1", 0, 3)
hh_2 = makeTH1(trnd, "hh2", -3, 1)
hh_3 = makeTH1(trnd, "hh3", +3, 4)

# Declare observable x
x = ROOT.RooRealVar("x", "x", -10, 10)
Expand All @@ -76,7 +77,7 @@ def makeTTree():
# Importing a ROOT TTree into a RooDataSet with cuts
# --------------------------------------------------------------------

tree = makeTTree()
tree = makeTTree(trnd)

# Define observables y,z
y = ROOT.RooRealVar("y", "y", -10, 10)
Expand Down
7 changes: 5 additions & 2 deletions tutorials/roofit/rf609_xychi2fit.C
Original file line number Diff line number Diff line change
Expand Up @@ -37,18 +37,21 @@ void rf609_xychi2fit()
RooRealVar y("y", "y", -10, 200);
RooDataSet dxy("dxy", "dxy", {x, y}, StoreError({x, y}));

TRandom3 trnd{};

// Fill an example dataset with X,err(X),Y,err(Y) values
RooArgSet vars{x, y};
for (int i = 0; i <= 10; i++) {

// Set X value and error
x = -10 + 2 * i;
x.setError(i < 5 ? 0.5 / 1. : 1.0 / 1.);

// Set Y value and error
y = x.getVal() * x.getVal() + 4 * fabs(gRandom->Gaus());
y = x.getVal() * x.getVal() + 4 * std::abs(trnd.Gaus());
y.setError(sqrt(y.getVal()));

dxy.add({x, y});
dxy.add(vars);
}

// P e r f o r m c h i 2 f i t t o X + / - d x a n d Y + / - d Y v a l u e s
Expand Down
3 changes: 2 additions & 1 deletion tutorials/roofit/rf609_xychi2fit.py
Original file line number Diff line number Diff line change
Expand Up @@ -14,6 +14,7 @@
import ROOT
import math

trnd = ROOT.TRandom3()

# Create dataset with X and Y values
# -------------------------------------------------------------------
Expand All @@ -35,7 +36,7 @@
x.setError((0.5 / 1.0) if (i < 5) else (1.0 / 1.0))

# Set Y value and error
y.setVal(x.getVal() * x.getVal() + 4 * abs(ROOT.gRandom.Gaus()))
y.setVal(x.getVal() * x.getVal() + 4 * abs(trnd.Gaus()))
y.setError(math.sqrt(y.getVal()))

dxy.add({x, y})
Expand Down

0 comments on commit ab46cb8

Please sign in to comment.