Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
rf408_RDataFrameToRooFit.py
Go to the documentation of this file.
1## \file
2## \ingroup tutorial_roofit
3## \notebook
4## Fill RooDataSet/RooDataHist in RDataFrame.
5##
6## This tutorial shows how to fill RooFit data classes directly from RDataFrame.
7## Using two small helpers, we tell RDataFrame where the data has to go.
8##
9## \macro_code
10## \macro_output
11##
12## \date July 2021
13## \author Harshal Shende, Stephan Hageboeck (C++ version)
14
15import ROOT
16import math
17
18
19# Set up
20# ------------------------
21
22# We create an RDataFrame with two columns filled with 2 million random numbers.
23df = ROOT.RDataFrame(2000000).Define("x", "gRandom->Uniform(-5., 5.)").Define("y", "gRandom->Gaus(1., 3.)")
24
25
26# We create RooFit variables that will represent the dataset.
27x = ROOT.RooRealVar("x", "x", -5.0, 5.0)
28y = ROOT.RooRealVar("y", "y", -50.0, 50.0)
29x.setBins(10)
30y.setBins(20)
31
32
33# Booking the creation of RooDataSet / RooDataHist in RDataFrame
34# ----------------------------------------------------------------
35
36# Method 1:
37# ---------
38# We directly book the RooDataSetHelper action.
39# We need to pass
40# - the RDataFrame column types as template parameters
41# - the constructor arguments for RooDataSet (they follow the same syntax as the usual RooDataSet constructors)
42# - the column names that RDataFrame should fill into the dataset
43
44# NOTE: RDataFrame columns are matched to RooFit variables by position, *not by name*!
45#
46# The returned object is not yet a RooDataSet, but an RResultPtr that will be
47# lazy-evaluated once you call GetValue() on it. We will only evaluate the
48# RResultPtr once all other RDataFrame related actions are declared. This way
49# we trigger the event loop computation only once, which will improve the
50# runtime significantly.
51#
52# To learn more about lazy actions, see:
53# https://root.cern/doc/master/classROOT_1_1RDataFrame.html#actions
54roo_data_set_result = df.Book(
55 ROOT.std.move(ROOT.RooDataSetHelper("dataset", "Title of dataset", ROOT.RooArgSet(x, y))), ("x", "y")
56)
57
58# Method 2:
59# ---------
60# We first declare the RooDataHistHelper
61rdhMaker = ROOT.RooDataHistHelper("dataset", "Title of dataset", ROOT.RooArgSet(x, y))
62
63# Then, we move it into an RDataFrame action:
64roo_data_hist_result = df.Book(ROOT.std.move(rdhMaker), ("x", "y"))
65
66
67# Run it and inspect the results
68# -------------------------------
69
70# At this point, all RDF actions were defined (namely, the `Book` operations),
71# so we can get values from the RResultPtr objects, triggering the event loop
72# and getting the actual RooFit data objects.
73roo_data_set = roo_data_set_result.GetValue()
74roo_data_hist = roo_data_hist_result.GetValue()
75
76# Let's inspect the dataset / datahist.
77
78
79def print_data(data):
80 print("")
81 data.Print()
82 for i in range(min(data.numEntries(), 20)):
83 print(
84 "("
85 + ", ".join(["{0:8.3f}".format(var.getVal()) for var in data.get(i)])
86 + ", ) weight={0:10.3f}".format(data.weight())
87 )
88
89 print("mean(x) = {0:.3f}".format(data.mean(x)) + "\tsigma(x) = {0:.3f}".format(math.sqrt(data.moment(x, 2.0))))
90 print("mean(y) = {0:.3f}".format(data.mean(y)) + "\tsigma(y) = {0:.3f}\n".format(math.sqrt(data.moment(y, 2.0))))
91
92
93print_data(roo_data_set)
94print_data(roo_data_hist)
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void char Point_t Rectangle_t WindowAttributes_t Float_t Float_t Float_t Int_t Int_t UInt_t UInt_t Rectangle_t Int_t Int_t Window_t TString Int_t GCValues_t GetPrimarySelectionOwner GetDisplay GetScreen GetColormap GetNativeEvent const char const char dpyName wid window const char font_name cursor keysym reg const char only_if_exist regb h Point_t winding char text const char depth char const char Int_t count const char ColorStruct_t color const char Pixmap_t Pixmap_t PictureAttributes_t attr const char char ret_data h unsigned char height h Atom_t Int_t ULong_t ULong_t unsigned char prop_list Atom_t Atom_t Atom_t Time_t format
ROOT's RDataFrame offers a modern, high-level interface for analysis of data stored in TTree ,...