rf401_importttreethx

'DATA AND CATEGORIES' RooFit tutorial macro #401

Overview of advanced option for importing data from ROOT ROOT.TTree and ROOT.THx histograms Basic import options are demonstrated in rf102_dataimport.py

Author: Clemens Lange, Wouter Verkerke (C version)
This notebook tutorial was automatically generated with ROOTBOOK-izer from the macro found in the ROOT repository on Sunday, November 27, 2022 at 11:07 AM.

In [1]:
import ROOT
from array import array


def makeTH1(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))

    return hh


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

    tree = ROOT.TTree("tree", "tree")
    px = array("d", [0])
    py = array("d", [0])
    pz = array("d", [0])
    pi = array("i", [0])
    tree.Branch("x", px, "x/D")
    tree.Branch("y", py, "y/D")
    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)
        pi[0] = i % 3
        tree.Fill()

    return tree
Welcome to JupyROOT 6.27/01

Import multiple TH1 into a RooDataHist

Create thee ROOT ROOT.TH1 histograms

In [2]:
hh_1 = makeTH1("hh1", 0, 3)
hh_2 = makeTH1("hh2", -3, 1)
hh_3 = makeTH1("hh3", +3, 4)

Declare observable x

In [3]:
x = ROOT.RooRealVar("x", "x", -10, 10)

Create category observable c that serves as index for the ROOT histograms

In [4]:
c = ROOT.RooCategory("c", "c")
c.defineType("SampleA")
c.defineType("SampleB")
c.defineType("SampleC")
Out[4]:
False

Create a binned dataset that imports contents of all ROOT.TH1 mapped by index category c

In [5]:
dh = ROOT.RooDataHist("dh", "dh", [x], Index=c, Import={"SampleA": hh_1, "SampleB": hh_2, "SampleC": hh_3})
dh.Print()

dh2 = ROOT.RooDataHist("dh", "dh", [x], Index=c, Import={"SampleA": hh_1, "SampleB": hh_2, "SampleC": hh_3})
dh2.Print()
RooDataHist::dh[c,x] = 300 bins (2964 weights)
RooDataHist::dh[c,x] = 300 bins (2964 weights)

Importing a ROOT TTree into a RooDataSet with cuts

In [6]:
tree = makeTTree()

Define observables y,z

In [7]:
y = ROOT.RooRealVar("y", "y", -10, 10)
z = ROOT.RooRealVar("z", "z", -10, 10)

Import only observables (y,z)

In [8]:
ds = ROOT.RooDataSet("ds", "ds", {x, y}, Import=tree)
ds.Print()
[#1] INFO:DataHandling -- RooTreeDataStore::loadValues(ds) Skipping event #7 because y cannot accommodate the value 13.3845
[#1] INFO:DataHandling -- RooTreeDataStore::loadValues(ds) Skipping event #8 because y cannot accommodate the value 11.1861
[#1] INFO:DataHandling -- RooTreeDataStore::loadValues(ds) Skipping event #12 because y cannot accommodate the value 13.7009
[#1] INFO:DataHandling -- RooTreeDataStore::loadValues(ds) Skipping event #14 because y cannot accommodate the value -10.6852
[#1] INFO:DataHandling -- RooTreeDataStore::loadValues(ds) Skipping ...
[#0] WARNING:DataHandling -- RooTreeDataStore::loadValues(ds) Ignored 35 out-of-range events
RooDataSet::ds[y,x] = 65 entries

Import observables (x,y,z) but only event for which (y+z<0) is ROOT.True Import observables (x,y,z) but only event for which (y+z<0) is ROOT.True

In [9]:
ds2 = ROOT.RooDataSet("ds2", "ds2", {x, y, z}, Import=tree, Cut="y+z<0")
ds2.Print()
[#1] INFO:InputArguments -- The formula y+z<0 claims to use the variables (y,x,z) but only (y,z) seem to be in use.
  inputs:         y+z<0
[#1] INFO:DataHandling -- RooTreeDataStore::loadValues(ds2) Skipping event #7 because y cannot accommodate the value 13.3845
[#1] INFO:DataHandling -- RooTreeDataStore::loadValues(ds2) Skipping event #8 because y cannot accommodate the value 11.1861
[#1] INFO:DataHandling -- RooTreeDataStore::loadValues(ds2) Skipping event #12 because y cannot accommodate the value 13.7009
[#1] INFO:DataHandling -- RooTreeDataStore::loadValues(ds2) Skipping event #14 because y cannot accommodate the value -10.6852
[#1] INFO:DataHandling -- RooTreeDataStore::loadValues(ds2) Skipping ...
[#0] WARNING:DataHandling -- RooTreeDataStore::loadValues(ds2) Ignored 36 out-of-range events
RooDataSet::ds2[y,x,z] = 26 entries

Importing integer ROOT TTree branches

Import integer tree branch as ROOT.RooRealVar

In [10]:
i = ROOT.RooRealVar("i", "i", 0, 5)
ds3 = ROOT.RooDataSet("ds3", "ds3", {i, x}, Import=tree)
ds3.Print()
[#1] INFO:DataHandling -- RooAbsReal::attachToTree(i) TTree Int_t branch i will be converted to double precision.
RooDataSet::ds3[i,x] = 100 entries

Define category i

In [11]:
icat = ROOT.RooCategory("i", "i", {"State0": 0, "State1": 1})

Import integer tree branch as ROOT.RooCategory (only events with i==0 and i==1 will be imported as those are the only defined states)

In [12]:
ds4 = ROOT.RooDataSet("ds4", "ds4", {icat, x}, Import=tree)
ds4.Print()
[#1] INFO:DataHandling -- RooTreeDataStore::loadValues(ds4) Skipping event #2 because i cannot accommodate the value 0
[#1] INFO:DataHandling -- RooTreeDataStore::loadValues(ds4) Skipping event #5 because i cannot accommodate the value 0
[#1] INFO:DataHandling -- RooTreeDataStore::loadValues(ds4) Skipping event #8 because i cannot accommodate the value 0
[#1] INFO:DataHandling -- RooTreeDataStore::loadValues(ds4) Skipping event #11 because i cannot accommodate the value 0
[#1] INFO:DataHandling -- RooTreeDataStore::loadValues(ds4) Skipping ...
[#0] WARNING:DataHandling -- RooTreeDataStore::loadValues(ds4) Ignored 33 out-of-range events
RooDataSet::ds4[i,x] = 67 entries

Import multiple RooDataSets into a RooDataSet

Create three ROOT.RooDataSets in (y,z)

In [13]:
dsA = ds2.reduce({x, y}, "z<-5")
dsB = ds2.reduce({x, y}, "abs(z)<5")
dsC = ds2.reduce({x, y}, "z>5")

Create a dataset that imports contents of all the above datasets mapped by index category c

In [14]:
dsABC = ROOT.RooDataSet("dsABC", "dsABC", {x, y}, Index=c, Import={"SampleA": dsA, "SampleB": dsB, "SampleC": dsC})

dsABC.Print()
RooDataSet::dsABC[c,y,x] = 26 entries

Draw all canvases

In [15]:
from ROOT import gROOT 
gROOT.GetListOfCanvases().Draw()