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