Logo ROOT  
Reference Guide
rf506_msgservice.py
Go to the documentation of this file.
1## \file
2## \ingroup tutorial_roofit
3## \notebook -nodraw
4##
5## Organization and simultaneous fits: tuning and customizing the ROOT.RooFit message logging facility
6##
7## \macro_code
8##
9## \date February 2018
10## \author Clemens Lange, Wouter Verkerke (C++ version)
11
12import ROOT
13
14# Create pdf
15# --------------------
16
17# Construct gauss(x,m,s)
18x = ROOT.RooRealVar("x", "x", -10, 10)
19m = ROOT.RooRealVar("m", "m", 0, -10, 10)
20s = ROOT.RooRealVar("s", "s", 1, -10, 10)
21gauss = ROOT.RooGaussian("g", "g", x, m, s)
22
23# Construct poly(x,p0)
24p0 = ROOT.RooRealVar("p0", "p0", 0.01, 0., 1.)
25poly = ROOT.RooPolynomial("p", "p", x, ROOT.RooArgList(p0))
26
27# model = f*gauss(x) + (1-f)*poly(x)
28f = ROOT.RooRealVar("f", "f", 0.5, 0., 1.)
29model = ROOT.RooAddPdf("model", "model", ROOT.RooArgList(
30 gauss, poly), ROOT.RooArgList(f))
31
32data = model.generate(ROOT.RooArgSet(x), 10)
33
34# Print configuration of message service
35# ------------------------------------------
36
37# Print streams configuration
38ROOT.RooMsgService.instance().Print()
39
40# Adding integration topic to existing INFO stream
41# ---------------------------------------------------
42
43# Print streams configuration
44ROOT.RooMsgService.instance().Print()
45
46# Add Integration topic to existing INFO stream
47ROOT.RooMsgService.instance().getStream(1).addTopic(ROOT.RooFit.Integration)
48
49# Construct integral over gauss to demonstrate message stream
50igauss = gauss.createIntegral(ROOT.RooArgSet(x))
51igauss.Print()
52
53# Print streams configuration in verbose, also shows inactive streams
54ROOT.RooMsgService.instance().Print()
55
56# Remove stream
57ROOT.RooMsgService.instance().getStream(1).removeTopic(ROOT.RooFit.Integration)
58
59# Examples of pdf value tracing
60# -----------------------------------------------------------------------
61
62# Show DEBUG level message on function tracing, ROOT.RooGaussian only
63ROOT.RooMsgService.instance().addStream(
64 ROOT.RooFit.DEBUG,
65 ROOT.RooFit.Topic(
66 ROOT.RooFit.Tracing),
67 ROOT.RooFit.ClassName("RooGaussian"))
68
69# Perform a fit to generate some tracing messages
70model.fitTo(data, ROOT.RooFit.Verbose(ROOT.kTRUE))
71
72# Reset message service to default stream configuration
73ROOT.RooMsgService.instance().reset()
74
75# Show DEBUG level message on function tracing on all objects, output to
76# file
77ROOT.RooMsgService.instance().addStream(
78 ROOT.RooFit.DEBUG,
79 ROOT.RooFit.Topic(
80 ROOT.RooFit.Tracing),
81 ROOT.RooFit.OutputFile("rf506_debug.log"))
82
83# Perform a fit to generate some tracing messages
84model.fitTo(data, ROOT.RooFit.Verbose(ROOT.kTRUE))
85
86# Reset message service to default stream configuration
87ROOT.RooMsgService.instance().reset()
88
89# Example of another debugging stream
90# ---------------------------------------------------------------------
91
92# Show DEBUG level messages on client/server link state management
93ROOT.RooMsgService.instance().addStream(
94 ROOT.RooFit.DEBUG, ROOT.RooFit.Topic(ROOT.RooFit.LinkStateMgmt))
95ROOT.RooMsgService.instance().Print("v")
96
97# Clone composite pdf g to trigger some link state management activity
98gprime = gauss.cloneTree()
99gprime.Print()
100
101# Reset message service to default stream configuration
102ROOT.RooMsgService.instance().reset()
void Print(std::ostream &os, const OptionType &opt)