Organisation and simultaneous fits: tuning and customizing the RooFit message logging facility
Author: Wouter Verkerke
This notebook tutorial was automatically generated with ROOTBOOK-izer from the macro found in the ROOT repository on Tuesday, March 19, 2024 at 07:16 PM.
%%cpp -d
#include "RooRealVar.h"
#include "RooDataSet.h"
#include "RooGaussian.h"
#include "RooPolynomial.h"
#include "RooAddPdf.h"
#include "TCanvas.h"
#include "TAxis.h"
#include "RooPlot.h"
#include "RooMsgService.h"
using namespace RooFit;
Construct gauss(x,m,s)
RooRealVar x("x", "x", -10, 10);
RooRealVar m("m", "m", 0, -10, 10);
RooRealVar s("s", "s", 1, -10, 10);
RooGaussian gauss("g", "g", x, m, s);
[#0] WARNING:InputArguments -- The parameter 's' with range [-10, 10] of the RooGaussian 'g' exceeds the safe range of (0, inf). Advise to limit its range.
Construct poly(x,p0)
RooRealVar p0("p0", "p0", 0.01, 0., 1.);
RooPolynomial poly("p", "p", x, p0);
Construct model = fgauss(x) + (1-f)poly(x)
RooRealVar f("f", "f", 0.5, 0., 1.);
RooAddPdf model("model", "model", RooArgSet(gauss, poly), f);
std::unique_ptr<RooDataSet> data{model.generate(x, 10)};
input_line_50:5:1: warning: 'data' shadows a declaration with the same name in the 'std' namespace; use '::data' to reference this declaration std::unique_ptr<RooDataSet> data{model.generate(x, 10)}; ^
Print streams configuration
RooMsgService::instance().Print();
cout << endl;
Active Message streams [0] MinLevel = PROGRESS Topic = Generation Minimization Plotting Fitting Integration LinkStateMgmt Eval Caching Optimization ObjectHandling InputArguments Tracing Contents DataHandling NumericIntegration FastEvaluations [1] MinLevel = INFO Topic = Minimization Plotting Fitting Eval Caching ObjectHandling InputArguments DataHandling NumericIntegration [2] MinLevel = INFO Topic = HistFactory
Print streams configuration
RooMsgService::instance().Print();
cout << endl;
Active Message streams [0] MinLevel = PROGRESS Topic = Generation Minimization Plotting Fitting Integration LinkStateMgmt Eval Caching Optimization ObjectHandling InputArguments Tracing Contents DataHandling NumericIntegration FastEvaluations [1] MinLevel = INFO Topic = Minimization Plotting Fitting Eval Caching ObjectHandling InputArguments DataHandling NumericIntegration [2] MinLevel = INFO Topic = HistFactory
Add Integration topic to existing INFO stream
RooMsgService::instance().getStream(1).addTopic(Integration);
Construct integral over gauss to demonstrate new message stream
std::unique_ptr<RooAbsReal> igauss{gauss.createIntegral(x)};
igauss->Print();
[#1] INFO:Integration -- RooRealIntegral::ctor(g_Int[x]) Constructing integral of function g over observables(x) with normalization () with range identifier <none> [#1] INFO:Integration -- g: Observable x is suitable for analytical integration (if supported by p.d.f) [#1] INFO:Integration -- g: Function integrated observables (x) internally with code 1 [#1] INFO:Integration -- g: Observables (x) are analytically integrated with code 1 RooRealIntegral::g_Int[x][ Int gd[Ana](x) ] = 2.50663
Print streams configuration in verbose, which also shows inactive streams
cout << endl;
RooMsgService::instance().Print();
cout << endl;
Active Message streams [0] MinLevel = PROGRESS Topic = Generation Minimization Plotting Fitting Integration LinkStateMgmt Eval Caching Optimization ObjectHandling InputArguments Tracing Contents DataHandling NumericIntegration FastEvaluations [1] MinLevel = INFO Topic = Minimization Plotting Fitting Integration Eval Caching ObjectHandling InputArguments DataHandling NumericIntegration [2] MinLevel = INFO Topic = HistFactory
Remove stream
RooMsgService::instance().getStream(1).removeTopic(Integration);
Show DEBUG level message on function tracing, trace RooGaussian only
RooMsgService::instance().addStream(DEBUG, Topic(Tracing), ClassName("RooGaussian"));
Perform a fit to generate some tracing messages
model.fitTo(*data, Verbose(true));
input_line_58:2:15: error: reference to 'data' is ambiguous model.fitTo(*data, Verbose(true)); ^ input_line_50:5:29: note: candidate found by name lookup is 'data' std::unique_ptr<RooDataSet> data{model.generate(x, 10)}; ^ /usr/include/c++/9/bits/range_access.h:318:5: note: candidate found by name lookup is 'std::data' data(initializer_list<_Tp> __il) noexcept ^ /usr/include/c++/9/bits/range_access.h:289:5: note: candidate found by name lookup is 'std::data' data(_Container& __cont) noexcept(noexcept(__cont.data())) ^ /usr/include/c++/9/bits/range_access.h:299:5: note: candidate found by name lookup is 'std::data' data(const _Container& __cont) noexcept(noexcept(__cont.data())) ^ /usr/include/c++/9/bits/range_access.h:309:5: note: candidate found by name lookup is 'std::data' data(_Tp (&__array)[_Nm]) noexcept ^
Reset message service to default stream configuration
RooMsgService::instance().reset();
Show DEBUG level message on function tracing on all objects, redirect output to file
RooMsgService::instance().addStream(DEBUG, Topic(Tracing), OutputFile("rf506_debug.log"));
Perform a fit to generate some tracing messages
model.fitTo(*data, Verbose(true));
input_line_61:2:15: error: reference to 'data' is ambiguous model.fitTo(*data, Verbose(true)); ^ input_line_50:5:29: note: candidate found by name lookup is 'data' std::unique_ptr<RooDataSet> data{model.generate(x, 10)}; ^ /usr/include/c++/9/bits/range_access.h:318:5: note: candidate found by name lookup is 'std::data' data(initializer_list<_Tp> __il) noexcept ^ /usr/include/c++/9/bits/range_access.h:289:5: note: candidate found by name lookup is 'std::data' data(_Container& __cont) noexcept(noexcept(__cont.data())) ^ /usr/include/c++/9/bits/range_access.h:299:5: note: candidate found by name lookup is 'std::data' data(const _Container& __cont) noexcept(noexcept(__cont.data())) ^ /usr/include/c++/9/bits/range_access.h:309:5: note: candidate found by name lookup is 'std::data' data(_Tp (&__array)[_Nm]) noexcept ^
Reset message service to default stream configuration
RooMsgService::instance().reset();
Show DEBUG level messages on client/server link state management
RooMsgService::instance().addStream(DEBUG, Topic(LinkStateMgmt));
RooMsgService::instance().Print("v");
All Message streams [0] MinLevel = PROGRESS Topic = Generation Minimization Plotting Fitting Integration LinkStateMgmt Eval Caching Optimization ObjectHandling InputArguments Tracing Contents DataHandling NumericIntegration FastEvaluations [1] MinLevel = INFO Topic = Minimization Plotting Fitting Eval Caching ObjectHandling InputArguments DataHandling NumericIntegration [2] MinLevel = INFO Topic = HistFactory [3] MinLevel = DEBUG Topic = LinkStateMgmt
Clone composite pdf g to trigger some link state management activity
RooAbsArg *gprime = gauss.cloneTree();
gprime->Print();
[#3] DEBUG:LinkStateMgmt -- RooAbsArg::addServer(0x7fdd25d46240,g): adding server x(0x7fdd40009000) for value [#3] DEBUG:LinkStateMgmt -- RooAbsArg::addServer(0x7fdd25d46240,g): adding server m(0x7fdd400093e8) for value [#3] DEBUG:LinkStateMgmt -- RooAbsArg::addServer(0x7fdd25d46240,g): adding server s(0x7fdd400097d0) for value [#3] DEBUG:LinkStateMgmt -- RooAbsArg::addServer(0x7fdd25d46240,g): adding server x(0x7fdd40009000) for value [#3] DEBUG:LinkStateMgmt -- RooAbsArg::addServer(0x7fdd25d46240,g): adding server m(0x7fdd400093e8) for value [#3] DEBUG:LinkStateMgmt -- RooAbsArg::addServer(0x7fdd25d46240,g): adding server s(0x7fdd400097d0) for value RooGaussian::g[ x=x mean=m sigma=s ] = 1
Reset message service to default stream configuration
RooMsgService::instance().reset();