Logo ROOT  
Reference Guide
rf404_categories.py
Go to the documentation of this file.
1## \file
2## \ingroup tutorial_roofit
3## \notebook -nodraw
4##
5## Data and categories: working with ROOT.RooCategory objects to describe discrete variables
6##
7## \macro_code
8##
9## \date February 2018
10## \author Clemens Lange, Wouter Verkerke (C++ version)
11
12from __future__ import print_function
13import ROOT
14
15
16# Construct a category with labels
17# --------------------------------------------
18
19# Define a category with labels only
20tagCat = ROOT.RooCategory("tagCat", "Tagging category")
21tagCat.defineType("Lepton")
22tagCat.defineType("Kaon")
23tagCat.defineType("NetTagger-1")
24tagCat.defineType("NetTagger-2")
25tagCat.Print()
26
27# Construct a category with labels and indices
28# ------------------------------------------------
29
30# Define a category with explicitly numbered states
31b0flav = ROOT.RooCategory("b0flav", "B0 flavour eigenstate")
32b0flav.defineType("B0", -1)
33b0flav.defineType("B0bar", 1)
34b0flav.Print()
35
36# Generate dummy data for tabulation demo
37# ------------------------------------------------
38
39# Generate a dummy dataset
40x = ROOT.RooRealVar("x", "x", 0, 10)
41data = ROOT.RooPolynomial("p", "p", x).generate(
42 ROOT.RooArgSet(x, b0flav, tagCat), 10000)
43
44# Print tables of category contents of datasets
45# --------------------------------------------------
46
47# Tables are equivalent of plots for categories
48btable = data.table(b0flav)
49btable.Print()
50btable.Print("v")
51
52# Create table for subset of events matching cut expression
53ttable = data.table(tagCat, "x>8.23")
54ttable.Print()
55ttable.Print("v")
56
57# Create table for all (tagCat x b0flav) state combinations
58bttable = data.table(ROOT.RooArgSet(tagCat, b0flav))
59bttable.Print("v")
60
61# Retrieve number of events from table
62# Number can be non-integer if source dataset has weighed events
63nb0 = btable.get("B0")
64print("Number of events with B0 flavor is ", nb0)
65
66# Retrieve fraction of events with "Lepton" tag
67fracLep = ttable.getFrac("Lepton")
68print("Fraction of events tagged with Lepton tag is ", fracLep)
69
70# Defining ranges for plotting, fitting on categories
71# ------------------------------------------------------------------------------------------------------
72
73# Define named range as comma separated list of labels
74tagCat.setRange("good", "Lepton,Kaon")
75
76# Or add state names one by one
77tagCat.addToRange("soso", "NetTagger-1")
78tagCat.addToRange("soso", "NetTagger-2")
79
80# Use category range in dataset reduction specification
81goodData = data.reduce(ROOT.RooFit.CutRange("good"))
82goodData.table(tagCat).Print("v")
void Print(std::ostream &os, const OptionType &opt)