Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
rf404_categories.py
Go to the documentation of this file.
1## \file
2## \ingroup tutorial_roofit
3## \notebook -nodraw
4## Data and categories: working with ROOT.RooCategory objects to describe discrete variables
5##
6## \macro_code
7## \macro_output
8##
9## \date February 2018
10## \authors 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", {"B0": -1, "B0bar": 1})
32b0flav.Print()
33
34# Generate dummy data for tabulation demo
35# ------------------------------------------------
36
37# Generate a dummy dataset
38x = ROOT.RooRealVar("x", "x", 0, 10)
39data = ROOT.RooPolynomial("p", "p", x).generate({x, b0flav, tagCat}, 10000)
40
41# Print tables of category contents of datasets
42# --------------------------------------------------
43
44# Tables are equivalent of plots for categories
45btable = data.table(b0flav)
46btable.Print()
47btable.Print("v")
48
49# Create table for subset of events matching cut expression
50ttable = data.table(tagCat, "x>8.23")
51ttable.Print()
52ttable.Print("v")
53
54# Create table for all (tagCat x b0flav) state combinations
55bttable = data.table({tagCat, b0flav})
56bttable.Print("v")
57
58# Retrieve number of events from table
59# Number can be non-integer if source dataset has weighed events
60nb0 = btable.get("B0")
61print("Number of events with B0 flavor is ", nb0)
62
63# Retrieve fraction of events with "Lepton" tag
64fracLep = ttable.getFrac("Lepton")
65print("Fraction of events tagged with Lepton tag is ", fracLep)
66
67# Defining ranges for plotting, fitting on categories
68# ------------------------------------------------------------------------------------------------------
69
70# Define named range as comma separated list of labels
71tagCat.setRange("good", "Lepton,Kaon")
72
73# Or add state names one by one
74tagCat.addToRange("soso", "NetTagger-1")
75tagCat.addToRange("soso", "NetTagger-2")
76
77# Use category range in dataset reduction specification
78goodData = data.reduce(CutRange="good")
79goodData.table(tagCat).Print("v")
void Print(GNN_Data &d, std::string txt="")