Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
RooMinimizer.cxx
Go to the documentation of this file.
1/*****************************************************************************
2 * Project: RooFit *
3 * Package: RooFitCore *
4 * @(#)root/roofitcore:$Id$
5 * Authors: *
6 * WV, Wouter Verkerke, UC Santa Barbara, verkerke@slac.stanford.edu *
7 * DK, David Kirkby, UC Irvine, dkirkby@uci.edu *
8 * AL, Alfio Lazzaro, INFN Milan, alfio.lazzaro@mi.infn.it *
9 * PB, Patrick Bos, NL eScience Center, p.bos@esciencecenter.nl *
10 * *
11 * Redistribution and use in source and binary forms, *
12 * with or without modification, are permitted according to the terms *
13 * listed in LICENSE (http://roofit.sourceforge.net/license.txt) *
14 *****************************************************************************/
15
16/**
17\file RooMinimizer.cxx
18\class RooMinimizer
19\ingroup Roofitcore
20
21RooMinimizer is a wrapper class around ROOT::Fit:Fitter that
22provides a seamless interface between the minimizer functionality
23and the native RooFit interface.
24By default the Minimizer is MINUIT for classic mode and MINUIT2
25for parallelized mode (activated with the `RooFit::ModularL(true)`
26parameter or by passing a RooRealL function as the minimization
27target).
28RooMinimizer can minimize any RooAbsReal function with respect to
29its parameters. Usual choices for minimization are RooNLLVar
30and RooChi2Var
31RooMinimizer has methods corresponding to MINUIT functions like
32hesse(), migrad(), minos() etc. In each of these function calls
33the state of the MINUIT engine is synchronized with the state
34of the RooFit variables: any change in variables, change
35in the constant status etc is forwarded to MINUIT prior to
36execution of the MINUIT call. Afterwards the RooFit objects
37are resynchronized with the output state of MINUIT: changes
38parameter values, errors are propagated.
39Various methods are available to control verbosity, profiling,
40automatic PDF optimization.
41**/
42
43#include "RooMinimizer.h"
44
45#include "RooAbsMinimizerFcn.h"
46#include "RooArgSet.h"
47#include "RooArgList.h"
48#include "RooAbsReal.h"
49#include "RooDataSet.h"
50#include "RooRealVar.h"
51#include "RooSentinel.h"
52#include "RooMsgService.h"
53#include "RooPlot.h"
54#include "RooMinimizerFcn.h"
55#include "RooFitResult.h"
59#ifdef R__HAS_ROOFIT_MULTIPROCESS
62#endif
63
64#include "TClass.h"
65#include "Math/Minimizer.h"
66#include "TMarker.h"
67#include "TGraph.h"
68#include "Fit/FitConfig.h"
69
70#include <fstream>
71#include <iostream>
72#include <stdexcept> // logic_error
73
74using namespace std;
75
77;
78
79std::unique_ptr<ROOT::Fit::Fitter> RooMinimizer::_theFitter = {};
80
81////////////////////////////////////////////////////////////////////////////////
82/// Cleanup method called by atexit handler installed by RooSentinel
83/// to delete all global heap objects when the program is terminated
84
86{
87 _theFitter.reset();
88}
89
90////////////////////////////////////////////////////////////////////////////////
91/// Construct MINUIT interface to given function. Function can be anything,
92/// but is typically a -log(likelihood) implemented by RooNLLVar or a chi^2
93/// (implemented by RooChi2Var). Other frequent use cases are a RooAddition
94/// of a RooNLLVar plus a penalty or constraint term. This class propagates
95/// all RooFit information (floating parameters, their values and errors)
96/// to MINUIT before each MINUIT call and propagates all MINUIT information
97/// back to the RooFit object at the end of each call (updated parameter
98/// values, their (asymmetric errors) etc. The default MINUIT error level
99/// for HESSE and MINOS error analysis is taken from the defaultErrorLevel()
100/// value of the input function.
101
102/// Constructor that accepts all configuration in struct with RooAbsReal likelihood
103RooMinimizer::RooMinimizer(RooAbsReal &function, Config const &cfg) : _cfg(cfg)
104{
106 auto nll_real = dynamic_cast<RooFit::TestStatistics::RooRealL *>(&function);
107 if (nll_real != nullptr) {
108 if (_cfg.parallelize != 0) { // new test statistic with multiprocessing library with
109 // parallel likelihood or parallel gradient
110#ifdef R__HAS_ROOFIT_MULTIPROCESS
112 // Note that this is necessary because there is currently no serial-mode LikelihoodGradientWrapper.
113 // We intend to repurpose RooGradMinimizerFcn to build such a LikelihoodGradientSerial class.
114 coutI(InputArguments) << "Modular likelihood detected and likelihood parallelization requested, "
115 << "also setting parallel gradient calculation mode." << std::endl;
117 }
118 // If _cfg.parallelize is larger than zero set the number of workers to that value. Otherwise do not do
119 // anything and let RooFit::MultiProcess handle the number of workers
120 if (_cfg.parallelize > 0)
123
124 _fcn = std::make_unique<RooFit::TestStatistics::MinuitFcnGrad>(
125 nll_real->getRooAbsL(), this, _theFitter->Config().ParamsSettings(),
127 static_cast<RooFit::TestStatistics::LikelihoodMode>(int(_cfg.enableParallelDescent))},
129#else
130 throw std::logic_error(
131 "Parallel minimization requested, but ROOT was not compiled with multiprocessing enabled, "
132 "please recompile with -Droofit_multiprocess=ON for parallel evaluation");
133#endif
134 } else { // modular test statistic non parallel
135 coutW(InputArguments)
136 << "Requested modular likelihood without gradient parallelization, some features such as offsetting "
137 << "may not work yet. Non-modular likelihoods are more reliable without parallelization." << std::endl;
138 // The RooRealL that is used in the case where the modular likelihood is being passed to a RooMinimizerFcn does
139 // not have offsetting implemented. Therefore, offsetting will not work in this case. Other features might also
140 // not work since the RooRealL was not intended for minimization. Further development is required to make the
141 // MinuitFcnGrad also handle serial gradient minimization. The MinuitFcnGrad accepts a RooAbsL and has
142 // offsetting implemented, thus omitting the need for RooRealL minimization altogether.
143 _fcn = std::make_unique<RooMinimizerFcn>(&function, this);
144 }
145 } else {
146 if (_cfg.parallelize != 0) { // Old test statistic with parallel likelihood or gradient
147 throw std::logic_error("In RooMinimizer constructor: Selected likelihood evaluation but a "
148 "non-modular likelihood was given. Please supply ModularL(true) as an "
149 "argument to createNLL for modular likelihoods to use likelihood "
150 "or gradient parallelization.");
151 }
152 _fcn = std::make_unique<RooMinimizerFcn>(&function, this);
153 }
154 initMinimizerFcnDependentPart(function.defaultErrorLevel());
155};
156
157/// Initialize the part of the minimizer that is independent of the function to be minimized
159{
162
163 _theFitter = std::make_unique<ROOT::Fit::Fitter>();
164 _theFitter->Config().SetMinimizer(_cfg.minimizerType.c_str());
165 setEps(1.0); // default tolerance
166}
167
168/// Initialize the part of the minimizer that is dependent on the function to be minimized
170{
171 // default max number of calls
172 _theFitter->Config().MinimizerOptions().SetMaxIterations(500 * _fcn->getNDim());
173 _theFitter->Config().MinimizerOptions().SetMaxFunctionCalls(500 * _fcn->getNDim());
174
175 // Shut up for now
176 setPrintLevel(-1);
177
178 // Use +0.5 for 1-sigma errors
179 setErrorLevel(defaultErrorLevel);
180
181 // Declare our parameters to MINUIT
182 _fcn->Synchronize(_theFitter->Config().ParamsSettings());
183
184 // Now set default verbosity
185 if (RooMsgService::instance().silentMode()) {
186 setPrintLevel(-1);
187 } else {
188 setPrintLevel(1);
189 }
190
191 // Set user defined and default _fcn config
193
194 // Likelihood holds information on offsetting in old style, so do not set here unless explicitly set by user
195 if (_cfg.offsetting != -1) {
197 }
198}
199
200////////////////////////////////////////////////////////////////////////////////
201/// Destructor
202
204
205////////////////////////////////////////////////////////////////////////////////
206/// Change MINUIT strategy to istrat. Accepted codes
207/// are 0,1,2 and represent MINUIT strategies for dealing
208/// most efficiently with fast FCNs (0), expensive FCNs (2)
209/// and 'intermediate' FCNs (1)
210
212{
213 _theFitter->Config().MinimizerOptions().SetStrategy(istrat);
214}
215
216////////////////////////////////////////////////////////////////////////////////
217/// Change maximum number of MINUIT iterations
218/// (RooMinimizer default 500 * #%parameters)
219
221{
222 _theFitter->Config().MinimizerOptions().SetMaxIterations(n);
223}
224
225////////////////////////////////////////////////////////////////////////////////
226/// Change maximum number of likelihood function calss from MINUIT
227/// (RooMinimizer default 500 * #%parameters)
228
230{
231 _theFitter->Config().MinimizerOptions().SetMaxFunctionCalls(n);
232}
233
234////////////////////////////////////////////////////////////////////////////////
235/// Set the level for MINUIT error analysis to the given
236/// value. This function overrides the default value
237/// that is taken in the RooMinimizer constructor from
238/// the defaultErrorLevel() method of the input function
239
241{
242 _theFitter->Config().MinimizerOptions().SetErrorDef(level);
243}
244
245////////////////////////////////////////////////////////////////////////////////
246/// Change MINUIT epsilon
247
248void RooMinimizer::setEps(double eps)
249{
250 _theFitter->Config().MinimizerOptions().SetTolerance(eps);
251}
252
253////////////////////////////////////////////////////////////////////////////////
254/// Enable internal likelihood offsetting for enhanced numeric precision
255
257{
258 _cfg.offsetting = flag;
259 _fcn->setOffsetting(_cfg.offsetting);
260}
261
262////////////////////////////////////////////////////////////////////////////////
263/// Choose the minimizer algorithm.
264///
265/// Passing an empty string selects the default minimizer type returned by
266/// ROOT::Math::MinimizerOptions::DefaultMinimizerType().
267
268void RooMinimizer::setMinimizerType(std::string const &type)
269{
271
272 if ((_cfg.parallelize != 0) && _cfg.minimizerType != "Minuit2") {
273 std::stringstream ss;
274 ss << "In RooMinimizer::setMinimizerType: only Minuit2 is supported when not using classic function mode!";
275 if (type.empty()) {
276 ss << "\nPlease set it as your default minimizer via "
277 "ROOT::Math::MinimizerOptions::SetDefaultMinimizer(\"Minuit2\").";
278 }
279 throw std::invalid_argument(ss.str());
280 }
281}
282
283////////////////////////////////////////////////////////////////////////////////
284/// Return underlying ROOT fitter object
285
287{
288 return _theFitter.get();
289}
290
291////////////////////////////////////////////////////////////////////////////////
292/// Return underlying ROOT fitter object
293
295{
296 return _theFitter.get();
297}
298
300{
301 return _fcn->fit(*_theFitter);
302}
303
304////////////////////////////////////////////////////////////////////////////////
305/// Minimise the function passed in the constructor.
306/// \param[in] type Type of fitter to use, e.g. "Minuit" "Minuit2". Passing an
307/// empty string will select the default minimizer type of the
308/// RooMinimizer, as returned by
309/// ROOT::Math::MinimizerOptions::DefaultMinimizerType().
310/// \attention This overrides the default fitter of this RooMinimizer.
311/// \param[in] alg Fit algorithm to use. (Optional)
312int RooMinimizer::minimize(const char *type, const char *alg)
313{
314 if (_cfg.timingAnalysis)
315#ifdef R__HAS_ROOFIT_MULTIPROCESS
317#else
318 throw std::logic_error("ProcessTimer, but ROOT was not compiled with multiprocessing enabled, "
319 "please recompile with -Droofit_multiprocess=ON for logging with the "
320 "ProcessTimer.");
321#endif
322 _fcn->Synchronize(_theFitter->Config().ParamsSettings());
323
325 _theFitter->Config().SetMinimizer(type, alg);
326
327 profileStart();
328 {
329 auto ctx = makeEvalErrorContext();
330
331 bool ret = fitFcn();
332 _status = ((ret) ? _theFitter->Result().Status() : -1);
333 }
334 profileStop();
335 _fcn->BackProp(_theFitter->Result());
336
337 saveStatus("MINIMIZE", _status);
338
339 return _status;
340}
341
342////////////////////////////////////////////////////////////////////////////////
343/// Execute MIGRAD. Changes in parameter values
344/// and calculated errors are automatically
345/// propagated back the RooRealVars representing
346/// the floating parameters in the MINUIT operation.
347
349{
350 _fcn->Synchronize(_theFitter->Config().ParamsSettings());
351 profileStart();
352 {
353 auto ctx = makeEvalErrorContext();
354
355 _theFitter->Config().SetMinimizer(_cfg.minimizerType.c_str(), "migrad");
356 bool ret = fitFcn();
357 _status = ((ret) ? _theFitter->Result().Status() : -1);
358 }
359 profileStop();
360 _fcn->BackProp(_theFitter->Result());
361
362 saveStatus("MIGRAD", _status);
363
364 return _status;
365}
366
367////////////////////////////////////////////////////////////////////////////////
368/// Execute HESSE. Changes in parameter values
369/// and calculated errors are automatically
370/// propagated back the RooRealVars representing
371/// the floating parameters in the MINUIT operation.
372
374{
375 if (_theFitter->GetMinimizer() == 0) {
376 coutW(Minimization) << "RooMinimizer::hesse: Error, run Migrad before Hesse!" << endl;
377 _status = -1;
378 } else {
379
380 _fcn->Synchronize(_theFitter->Config().ParamsSettings());
381 profileStart();
382 {
383 auto ctx = makeEvalErrorContext();
384
385 _theFitter->Config().SetMinimizer(_cfg.minimizerType.c_str());
386 bool ret = _theFitter->CalculateHessErrors();
387 _status = ((ret) ? _theFitter->Result().Status() : -1);
388 }
389 profileStop();
390 _fcn->BackProp(_theFitter->Result());
391
392 saveStatus("HESSE", _status);
393 }
394
395 return _status;
396}
397
398////////////////////////////////////////////////////////////////////////////////
399/// Execute MINOS. Changes in parameter values
400/// and calculated errors are automatically
401/// propagated back the RooRealVars representing
402/// the floating parameters in the MINUIT operation.
403
405{
406 if (_theFitter->GetMinimizer() == 0) {
407 coutW(Minimization) << "RooMinimizer::minos: Error, run Migrad before Minos!" << endl;
408 _status = -1;
409 } else {
410
411 _fcn->Synchronize(_theFitter->Config().ParamsSettings());
412 profileStart();
413 {
414 auto ctx = makeEvalErrorContext();
415
416 _theFitter->Config().SetMinimizer(_cfg.minimizerType.c_str());
417 bool ret = _theFitter->CalculateMinosErrors();
418 _status = ((ret) ? _theFitter->Result().Status() : -1);
419 }
420
421 profileStop();
422 _fcn->BackProp(_theFitter->Result());
423
424 saveStatus("MINOS", _status);
425 }
426
427 return _status;
428}
429
430////////////////////////////////////////////////////////////////////////////////
431/// Execute MINOS for given list of parameters. Changes in parameter values
432/// and calculated errors are automatically
433/// propagated back the RooRealVars representing
434/// the floating parameters in the MINUIT operation.
435
436int RooMinimizer::minos(const RooArgSet &minosParamList)
437{
438 if (_theFitter->GetMinimizer() == 0) {
439 coutW(Minimization) << "RooMinimizer::minos: Error, run Migrad before Minos!" << endl;
440 _status = -1;
441 } else if (!minosParamList.empty()) {
442
443 _fcn->Synchronize(_theFitter->Config().ParamsSettings());
444 profileStart();
445 {
446 auto ctx = makeEvalErrorContext();
447
448 // get list of parameters for Minos
449 std::vector<unsigned int> paramInd;
450 for (RooAbsArg *arg : minosParamList) {
451 RooAbsArg *par = _fcn->GetFloatParamList()->find(arg->GetName());
452 if (par && !par->isConstant()) {
453 int index = _fcn->GetFloatParamList()->index(par);
454 paramInd.push_back(index);
455 }
456 }
457
458 if (paramInd.size()) {
459 // set the parameter indeces
460 _theFitter->Config().SetMinosErrors(paramInd);
461
462 _theFitter->Config().SetMinimizer(_cfg.minimizerType.c_str());
463 bool ret = _theFitter->CalculateMinosErrors();
464 _status = ((ret) ? _theFitter->Result().Status() : -1);
465 // to avoid that following minimization computes automatically the Minos errors
466 _theFitter->Config().SetMinosErrors(false);
467 }
468 }
469 profileStop();
470 _fcn->BackProp(_theFitter->Result());
471
472 saveStatus("MINOS", _status);
473 }
474
475 return _status;
476}
477
478////////////////////////////////////////////////////////////////////////////////
479/// Execute SEEK. Changes in parameter values
480/// and calculated errors are automatically
481/// propagated back the RooRealVars representing
482/// the floating parameters in the MINUIT operation.
483
485{
486 _fcn->Synchronize(_theFitter->Config().ParamsSettings());
487 profileStart();
488 {
489 auto ctx = makeEvalErrorContext();
490
491 _theFitter->Config().SetMinimizer(_cfg.minimizerType.c_str(), "seek");
492 bool ret = fitFcn();
493 _status = ((ret) ? _theFitter->Result().Status() : -1);
494 }
495 profileStop();
496 _fcn->BackProp(_theFitter->Result());
497
498 saveStatus("SEEK", _status);
499
500 return _status;
501}
502
503////////////////////////////////////////////////////////////////////////////////
504/// Execute SIMPLEX. Changes in parameter values
505/// and calculated errors are automatically
506/// propagated back the RooRealVars representing
507/// the floating parameters in the MINUIT operation.
508
510{
511 _fcn->Synchronize(_theFitter->Config().ParamsSettings());
512 profileStart();
513 {
514 auto ctx = makeEvalErrorContext();
515
516 _theFitter->Config().SetMinimizer(_cfg.minimizerType.c_str(), "simplex");
517 bool ret = fitFcn();
518 _status = ((ret) ? _theFitter->Result().Status() : -1);
519 }
520 profileStop();
521 _fcn->BackProp(_theFitter->Result());
522
523 saveStatus("SEEK", _status);
524
525 return _status;
526}
527
528////////////////////////////////////////////////////////////////////////////////
529/// Execute IMPROVE. Changes in parameter values
530/// and calculated errors are automatically
531/// propagated back the RooRealVars representing
532/// the floating parameters in the MINUIT operation.
533
535{
536 _fcn->Synchronize(_theFitter->Config().ParamsSettings());
537 profileStart();
538 {
539 auto ctx = makeEvalErrorContext();
540
541 _theFitter->Config().SetMinimizer(_cfg.minimizerType.c_str(), "migradimproved");
542 bool ret = fitFcn();
543 _status = ((ret) ? _theFitter->Result().Status() : -1);
544 }
545 profileStop();
546 _fcn->BackProp(_theFitter->Result());
547
548 saveStatus("IMPROVE", _status);
549
550 return _status;
551}
552
553////////////////////////////////////////////////////////////////////////////////
554/// Change the MINUIT internal printing level
555
557{
558 _theFitter->Config().MinimizerOptions().SetPrintLevel(newLevel + 1);
559}
560
561////////////////////////////////////////////////////////////////////////////////
562/// Get the MINUIT internal printing level
563
565{
566 return _theFitter->Config().MinimizerOptions().PrintLevel() + 1;
567}
568
569////////////////////////////////////////////////////////////////////////////////
570/// If flag is true, perform constant term optimization on
571/// function being minimized.
572
574{
575 _fcn->setOptimizeConst(flag);
576}
577
578////////////////////////////////////////////////////////////////////////////////
579/// Save and return a RooFitResult snapshot of current minimizer status.
580/// This snapshot contains the values of all constant parameters,
581/// the value of all floating parameters at RooMinimizer construction and
582/// after the last MINUIT operation, the MINUIT status, variance quality,
583/// EDM setting, number of calls with evaluation problems, the minimized
584/// function value and the full correlation matrix.
585
586RooFit::OwningPtr<RooFitResult> RooMinimizer::save(const char *userName, const char *userTitle)
587{
588 if (_theFitter->GetMinimizer() == 0) {
589 coutW(Minimization) << "RooMinimizer::save: Error, run minimization before!" << endl;
590 return nullptr;
591 }
592
593 TString name, title;
594 name = userName ? userName : Form("%s", _fcn->getFunctionName().c_str());
595 title = userTitle ? userTitle : Form("%s", _fcn->getFunctionTitle().c_str());
596 auto fitRes = std::make_unique<RooFitResult>(name, title);
597
598 // Move eventual fixed parameters in floatList to constList
599 RooArgList saveConstList(*(_fcn->GetConstParamList()));
600 RooArgList saveFloatInitList(*(_fcn->GetInitFloatParamList()));
601 RooArgList saveFloatFinalList(*(_fcn->GetFloatParamList()));
602 for (std::size_t i = 0; i < _fcn->GetFloatParamList()->size(); i++) {
603 RooAbsArg *par = _fcn->GetFloatParamList()->at(i);
604 if (par->isConstant()) {
605 saveFloatInitList.remove(*saveFloatInitList.find(par->GetName()), true);
606 saveFloatFinalList.remove(*par);
607 saveConstList.add(*par);
608 }
609 }
610 saveConstList.sort();
611
612 fitRes->setConstParList(saveConstList);
613 fitRes->setInitParList(saveFloatInitList);
614
615 // The fitter often clones the function. We therefore have to ask it for its copy.
616 const auto fitFcn = dynamic_cast<const RooAbsMinimizerFcn *>(_theFitter->GetFCN());
617 double removeOffset = 0.;
618 if (fitFcn) {
619 fitRes->setNumInvalidNLL(fitFcn->GetNumInvalidNLL());
620 removeOffset = -fitFcn->getOffset();
621 }
622
623 fitRes->setStatus(_status);
624 fitRes->setCovQual(_theFitter->GetMinimizer()->CovMatrixStatus());
625 fitRes->setMinNLL(_theFitter->Result().MinFcnValue() + removeOffset);
626 fitRes->setEDM(_theFitter->Result().Edm());
627 fitRes->setFinalParList(saveFloatFinalList);
628 if (!_extV) {
629 std::vector<double> globalCC;
630 TMatrixDSym corrs(_theFitter->Result().Parameters().size());
631 TMatrixDSym covs(_theFitter->Result().Parameters().size());
632 for (std::size_t ic = 0; ic < _theFitter->Result().Parameters().size(); ic++) {
633 globalCC.push_back(_theFitter->Result().GlobalCC(ic));
634 for (std::size_t ii = 0; ii < _theFitter->Result().Parameters().size(); ii++) {
635 corrs(ic, ii) = _theFitter->Result().Correlation(ic, ii);
636 covs(ic, ii) = _theFitter->Result().CovMatrix(ic, ii);
637 }
638 }
639 fitRes->fillCorrMatrix(globalCC, corrs, covs);
640 } else {
641 fitRes->setCovarianceMatrix(*_extV);
642 }
643
644 fitRes->setStatusHistory(_statusHistory);
645
646 return RooFit::Detail::owningPtr(std::move(fitRes));
647}
648
649////////////////////////////////////////////////////////////////////////////////
650/// Create and draw a TH2 with the error contours in the parameters `var1` and `var2`.
651/// \param[in] var1 The first parameter (x axis).
652/// \param[in] var2 The second parameter (y axis).
653/// \param[in] n1 First contour.
654/// \param[in] n2 Optional contour. 0 means don't draw.
655/// \param[in] n3 Optional contour. 0 means don't draw.
656/// \param[in] n4 Optional contour. 0 means don't draw.
657/// \param[in] n5 Optional contour. 0 means don't draw.
658/// \param[in] n6 Optional contour. 0 means don't draw.
659/// \param[in] npoints Number of points for evaluating the contour.
660///
661/// Up to six contours can be drawn using the arguments `n1` to `n6` to request the desired
662/// coverage in units of \f$ \sigma = n^2 \cdot \mathrm{ErrorDef} \f$.
663/// See ROOT::Math::Minimizer::ErrorDef().
664
665RooPlot *RooMinimizer::contour(RooRealVar &var1, RooRealVar &var2, double n1, double n2, double n3, double n4,
666 double n5, double n6, unsigned int npoints)
667{
668 RooArgList *params = _fcn->GetFloatParamList();
669 std::unique_ptr<RooArgList> paramSave{static_cast<RooArgList *>(params->snapshot())};
670
671 // Verify that both variables are floating parameters of PDF
672 int index1 = _fcn->GetFloatParamList()->index(&var1);
673 if (index1 < 0) {
674 coutE(Minimization) << "RooMinimizer::contour(" << GetName() << ") ERROR: " << var1.GetName()
675 << " is not a floating parameter of " << _fcn->getFunctionName() << endl;
676 return nullptr;
677 }
678
679 int index2 = _fcn->GetFloatParamList()->index(&var2);
680 if (index2 < 0) {
681 coutE(Minimization) << "RooMinimizer::contour(" << GetName() << ") ERROR: " << var2.GetName()
682 << " is not a floating parameter of PDF " << _fcn->getFunctionName() << endl;
683 return nullptr;
684 }
685
686 // create and draw a frame
687 RooPlot *frame = new RooPlot(var1, var2);
688
689 // draw a point at the current parameter values
690 TMarker *point = new TMarker(var1.getVal(), var2.getVal(), 8);
691 frame->addObject(point);
692
693 // check first if a inimizer is available. If not means
694 // the minimization is not done , so do it
695 if (_theFitter->GetMinimizer() == 0) {
696 coutW(Minimization) << "RooMinimizer::contour: Error, run Migrad before contours!" << endl;
697 return frame;
698 }
699
700 // remember our original value of ERRDEF
701 double errdef = _theFitter->GetMinimizer()->ErrorDef();
702
703 double n[6];
704 n[0] = n1;
705 n[1] = n2;
706 n[2] = n3;
707 n[3] = n4;
708 n[4] = n5;
709 n[5] = n6;
710
711 for (int ic = 0; ic < 6; ic++) {
712 if (n[ic] > 0) {
713
714 // set the value corresponding to an n1-sigma contour
715 _theFitter->GetMinimizer()->SetErrorDef(n[ic] * n[ic] * errdef);
716
717 // calculate and draw the contour
718 std::vector<double> xcoor(npoints + 1);
719 std::vector<double> ycoor(npoints + 1);
720 bool ret = _theFitter->GetMinimizer()->Contour(index1, index2, npoints, xcoor.data(), ycoor.data());
721
722 if (!ret) {
723 coutE(Minimization) << "RooMinimizer::contour(" << GetName()
724 << ") ERROR: MINUIT did not return a contour graph for n=" << n[ic] << endl;
725 } else {
726 xcoor[npoints] = xcoor[0];
727 ycoor[npoints] = ycoor[0];
728 TGraph *graph = new TGraph(npoints + 1, xcoor.data(), ycoor.data());
729
730 graph->SetName(Form("contour_%s_n%f", _fcn->getFunctionName().c_str(), n[ic]));
731 graph->SetLineStyle(ic + 1);
732 graph->SetLineWidth(2);
733 graph->SetLineColor(kBlue);
734 frame->addObject(graph, "L");
735 }
736 }
737 }
738
739 // restore the original ERRDEF
740 _theFitter->GetMinimizer()->SetErrorDef(errdef);
741
742 // restore parameter values
743 params->assign(*paramSave);
744
745 return frame;
746}
747
748////////////////////////////////////////////////////////////////////////////////
749/// Add parameters in metadata field to process timer
750
752{
753#ifdef R__HAS_ROOFIT_MULTIPROCESS
754 // parameter indices for use in timing heat matrix
755 std::vector<std::string> parameter_names;
756 for (auto &&parameter : *_fcn->GetFloatParamList()) {
757 parameter_names.push_back(parameter->GetName());
758 if (_cfg.verbose) {
759 coutI(Minimization) << "parameter name: " << parameter_names.back() << std::endl;
760 }
761 }
763#else
764 coutI(Minimization) << "Not adding parameters to processtimer because multiprocessing "
765 << "is not enabled." << std::endl;
766#endif
767}
768
769////////////////////////////////////////////////////////////////////////////////
770/// Start profiling timer
771
773{
774 if (_cfg.profile) {
775 _timer.Start();
776 _cumulTimer.Start(_profileStart ? false : true);
777 _profileStart = true;
778 }
779}
780
781////////////////////////////////////////////////////////////////////////////////
782/// Stop profiling timer and report results of last session
783
785{
786 if (_cfg.profile) {
787 _timer.Stop();
789 coutI(Minimization) << "Command timer: ";
790 _timer.Print();
791 coutI(Minimization) << "Session timer: ";
793 }
794}
795
797{
798 auto *fitterFcn = fitter()->GetFCN();
799 return fitterFcn ? fitterFcn : _fcn->getMultiGenFcn();
800}
801
802////////////////////////////////////////////////////////////////////////////////
803/// Apply results of given external covariance matrix. i.e. propagate its errors
804/// to all RRV parameter representations and give this matrix instead of the
805/// HESSE matrix at the next save() call
806
808{
809 _extV.reset(static_cast<TMatrixDSym *>(V.Clone()));
810 _fcn->ApplyCovarianceMatrix(*_extV);
811}
812
814{
816}
817
819{
820 // Import the results of the last fit performed, interpreting
821 // the fit parameters as the given varList of parameters.
822
823 if (_theFitter == 0 || _theFitter->GetMinimizer() == 0) {
824 oocoutE(nullptr, InputArguments) << "RooMinimizer::save: Error, run minimization before!" << endl;
825 return nullptr;
826 }
827
828 // Verify length of supplied varList
829 if (!varList.empty() && varList.size() != _theFitter->Result().NTotalParameters()) {
830 oocoutE(nullptr, InputArguments)
831 << "RooMinimizer::lastMinuitFit: ERROR: supplied variable list must be either empty " << endl
832 << " or match the number of variables of the last fit ("
833 << _theFitter->Result().NTotalParameters() << ")" << endl;
834 return nullptr;
835 }
836
837 // Verify that all members of varList are of type RooRealVar
838 for (RooAbsArg *arg : varList) {
839 if (!dynamic_cast<RooRealVar *>(arg)) {
840 oocoutE(nullptr, InputArguments) << "RooMinimizer::lastMinuitFit: ERROR: variable '" << arg->GetName()
841 << "' is not of type RooRealVar" << endl;
842 return nullptr;
843 }
844 }
845
846 auto res = std::make_unique<RooFitResult>("lastMinuitFit", "Last MINUIT fit");
847
848 // Extract names of fit parameters
849 // and construct corresponding RooRealVars
850 RooArgList constPars("constPars");
851 RooArgList floatPars("floatPars");
852
853 unsigned int i;
854 for (i = 0; i < _theFitter->Result().NTotalParameters(); ++i) {
855
856 TString varName(_theFitter->Result().GetParameterName(i));
857 bool isConst(_theFitter->Result().IsParameterFixed(i));
858
859 double xlo = _theFitter->Config().ParSettings(i).LowerLimit();
860 double xhi = _theFitter->Config().ParSettings(i).UpperLimit();
861 double xerr = _theFitter->Result().Error(i);
862 double xval = _theFitter->Result().Value(i);
863
864 std::unique_ptr<RooRealVar> var;
865 if (varList.empty()) {
866
867 if ((xlo < xhi) && !isConst) {
868 var = std::make_unique<RooRealVar>(varName, varName, xval, xlo, xhi);
869 } else {
870 var = std::make_unique<RooRealVar>(varName, varName, xval);
871 }
872 var->setConstant(isConst);
873 } else {
874
875 var.reset(static_cast<RooRealVar *>(varList.at(i)->Clone()));
876 var->setConstant(isConst);
877 var->setVal(xval);
878 if (xlo < xhi) {
879 var->setRange(xlo, xhi);
880 }
881
882 if (varName.CompareTo(var->GetName())) {
883 oocoutI(nullptr, Eval) << "RooMinimizer::lastMinuitFit: fit parameter '" << varName
884 << "' stored in variable '" << var->GetName() << "'" << endl;
885 }
886 }
887
888 if (isConst) {
889 constPars.addOwned(std::move(var));
890 } else {
891 var->setError(xerr);
892 floatPars.addOwned(std::move(var));
893 }
894 }
895
896 res->setConstParList(constPars);
897 res->setInitParList(floatPars);
898 res->setFinalParList(floatPars);
899 res->setMinNLL(_theFitter->Result().MinFcnValue());
900 res->setEDM(_theFitter->Result().Edm());
901 res->setCovQual(_theFitter->GetMinimizer()->CovMatrixStatus());
902 res->setStatus(_theFitter->Result().Status());
903 std::vector<double> globalCC;
904 TMatrixDSym corrs(_theFitter->Result().Parameters().size());
905 TMatrixDSym covs(_theFitter->Result().Parameters().size());
906 for (unsigned int ic = 0; ic < _theFitter->Result().Parameters().size(); ic++) {
907 globalCC.push_back(_theFitter->Result().GlobalCC(ic));
908 for (unsigned int ii = 0; ii < _theFitter->Result().Parameters().size(); ii++) {
909 corrs(ic, ii) = _theFitter->Result().Correlation(ic, ii);
910 covs(ic, ii) = _theFitter->Result().CovMatrix(ic, ii);
911 }
912 }
913 res->fillCorrMatrix(globalCC, corrs, covs);
914
915 return RooFit::Detail::owningPtr(std::move(res));
916}
917
918/// Try to recover from invalid function values. When invalid function values
919/// are encountered, a penalty term is returned to the minimiser to make it
920/// back off. This sets the strength of this penalty. \note A strength of zero
921/// is equivalent to a constant penalty (= the gradient vanishes, ROOT < 6.24).
922/// Positive values lead to a gradient pointing away from the undefined
923/// regions. Use ~10 to force the minimiser away from invalid function values.
925{
926 _cfg.recoverFromNaN = strength;
927}
928
929bool RooMinimizer::setLogFile(const char *logf)
930{
931 _cfg.logf = logf;
932 if (_cfg.logf)
933 return _fcn->SetLogFile(_cfg.logf);
934 else
935 return false;
936}
937
939{
940 return _fcn->evalCounter();
941}
943{
944 _fcn->zeroEvalCount();
945}
946
948{
949 return _fcn->getNDim();
950}
951
952std::ofstream *RooMinimizer::logfile()
953{
954 return _fcn->GetLogFile();
955}
957{
958 return _fcn->GetMaxFCN();
959}
960
962{
963#ifdef R__HAS_ROOFIT_MULTIPROCESS
965#else
966 return 0;
967#endif
968}
969
970std::unique_ptr<RooAbsReal::EvalErrorContext> RooMinimizer::makeEvalErrorContext() const
971{
973 // If evaluation error printing is disabled, we don't need to collect the
974 // errors and only need to count them. This significantly reduces the
975 // performance overhead when having evaluation errors.
977 return std::make_unique<RooAbsReal::EvalErrorContext>(m);
978}
size_t size(const MatrixT &matrix)
retrieve the size of a square matrix
#define coutI(a)
#define coutW(a)
#define oocoutE(o, a)
#define oocoutI(o, a)
#define coutE(a)
#define ClassImp(name)
Definition Rtypes.h:377
@ kBlue
Definition Rtypes.h:66
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void char Point_t Rectangle_t WindowAttributes_t index
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 type
char name[80]
Definition TGX11.cxx:110
char * Form(const char *fmt,...)
Formats a string in a circular formatting buffer.
Definition TString.cxx:2467
Fitter class, entry point for performing all type of fits.
Definition Fitter.h:77
ROOT::Math::IMultiGenFunction * GetFCN() const
return pointer to last used objective function (is NULL in case fit is not yet done) This pointer wil...
Definition Fitter.h:479
Documentation for the abstract class IBaseFunctionMultiDim.
Definition IFunction.h:62
static const std::string & DefaultMinimizerType()
RooAbsArg is the common abstract base class for objects that represent a value and a "shape" in RooFi...
Definition RooAbsArg.h:74
bool isConstant() const
Check if the "Constant" attribute is set.
Definition RooAbsArg.h:359
TObject * Clone(const char *newname=nullptr) const override
Make a clone of an object using the Streamer facility.
Definition RooAbsArg.h:86
virtual bool remove(const RooAbsArg &var, bool silent=false, bool matchByNameOnly=false)
Remove the specified argument from our list.
RooAbsCollection * snapshot(bool deepCopy=true) const
Take a snap shot of current collection contents.
virtual bool add(const RooAbsArg &var, bool silent=false)
Add the specified argument to list.
Int_t index(const RooAbsArg *arg) const
Returns index of given arg, or -1 if arg is not in the collection.
void assign(const RooAbsCollection &other) const
Sets the value, cache and constant attribute of any argument in our set that also appears in the othe...
Storage_t::size_type size() const
virtual bool addOwned(RooAbsArg &var, bool silent=false)
Add an argument and transfer the ownership to the collection.
void sort(bool reverse=false)
Sort collection using std::sort and name comparison.
RooAbsArg * find(const char *name) const
Find object with given name in list.
RooAbsReal is the common abstract base class for objects that represent a real value and implements f...
Definition RooAbsReal.h:62
double getVal(const RooArgSet *normalisationSet=nullptr) const
Evaluate object.
Definition RooAbsReal.h:91
static void clearEvalErrorLog()
Clear the stack of evaluation error messages.
RooArgList is a container object that can hold multiple RooAbsArg objects.
Definition RooArgList.h:22
RooAbsArg * at(Int_t idx) const
Return object at given index, or nullptr if index is out of range.
Definition RooArgList.h:110
RooArgSet is a container object that can hold multiple RooAbsArg objects.
Definition RooArgSet.h:55
static void setDefaultNWorkers(unsigned int N_workers)
Definition Config.cxx:67
static unsigned int getDefaultNWorkers()
Definition Config.cxx:92
static void setTimingAnalysis(bool timingAnalysis)
Definition Config.cxx:78
static void add_metadata(json data)
RooAbsReal that wraps RooAbsL likelihoods for use in RooFit outside of the RooMinimizer context.
Definition RooRealL.h:28
RooMinimizer is a wrapper class around ROOT::Fit:Fitter that provides a seamless interface between th...
void setRecoverFromNaNStrength(double strength)
Try to recover from invalid function values.
static int getPrintLevel()
Get the MINUIT internal printing level.
void optimizeConst(int flag)
If flag is true, perform constant term optimization on function being minimized.
void initMinimizerFirstPart()
Initialize the part of the minimizer that is independent of the function to be minimized.
std::ofstream * logfile()
int simplex()
Execute SIMPLEX.
std::unique_ptr< TMatrixDSym > _extV
void setMinimizerType(std::string const &type)
Choose the minimizer algorithm.
RooFit::OwningPtr< RooFitResult > save(const char *name=nullptr, const char *title=nullptr)
Save and return a RooFitResult snapshot of current minimizer status.
std::vector< std::pair< std::string, int > > _statusHistory
void profileStart()
Start profiling timer.
RooPlot * contour(RooRealVar &var1, RooRealVar &var2, double n1=1.0, double n2=2.0, double n3=0.0, double n4=0.0, double n5=0.0, double n6=0.0, unsigned int npoints=50)
Create and draw a TH2 with the error contours in the parameters var1 and var2.
bool setLogFile(const char *logf=nullptr)
void initMinimizerFcnDependentPart(double defaultErrorLevel)
Initialize the part of the minimizer that is dependent on the function to be minimized.
void profileStop()
Stop profiling timer and report results of last session.
int minos()
Execute MINOS.
double & maxFCN()
int hesse()
Execute HESSE.
void setErrorLevel(double level)
Set the level for MINUIT error analysis to the given value.
static std::unique_ptr< ROOT::Fit::Fitter > _theFitter
int migrad()
Execute MIGRAD.
int seek()
Execute SEEK.
void setEps(double eps)
Change MINUIT epsilon.
void setPrintLevel(int newLevel)
Change the MINUIT internal printing level.
int improve()
Execute IMPROVE.
static void cleanup()
Cleanup method called by atexit handler installed by RooSentinel to delete all global heap objects wh...
void setOffsetting(bool flag)
Enable internal likelihood offsetting for enhanced numeric precision.
TStopwatch _timer
RooMinimizer::Config _cfg
static RooFit::OwningPtr< RooFitResult > lastMinuitFit()
bool fitFcn() const
void saveStatus(const char *label, int status)
~RooMinimizer() override
Destructor.
int minimize(const char *type, const char *alg=nullptr)
Minimise the function passed in the constructor.
ROOT::Math::IMultiGenFunction * getMultiGenFcn() const
void setStrategy(int strat)
Change MINUIT strategy to istrat.
std::unique_ptr< RooAbsReal::EvalErrorContext > makeEvalErrorContext() const
RooMinimizer(RooAbsReal &function, Config const &cfg={})
Construct MINUIT interface to given function.
ROOT::Fit::Fitter * fitter()
Return underlying ROOT fitter object.
void setMaxFunctionCalls(int n)
Change maximum number of likelihood function calss from MINUIT (RooMinimizer default 500 * #parameter...
int evalCounter() const
TStopwatch _cumulTimer
int getNPar() const
void setMaxIterations(int n)
Change maximum number of MINUIT iterations (RooMinimizer default 500 * #parameters)
void addParamsToProcessTimer()
Add parameters in metadata field to process timer.
std::unique_ptr< RooAbsMinimizerFcn > _fcn
void applyCovarianceMatrix(TMatrixDSym const &V)
Apply results of given external covariance matrix.
static RooMsgService & instance()
Return reference to singleton instance.
A RooPlot is a plot frame and a container for graphics objects within that frame.
Definition RooPlot.h:43
void addObject(TObject *obj, Option_t *drawOptions="", bool invisible=false)
Add a generic object to this plot.
Definition RooPlot.cxx:380
RooRealVar represents a variable that can be changed from the outside.
Definition RooRealVar.h:40
static void activate()
Install atexit handler that calls CleanupRooFitAtExit() on program termination.
A TGraph is an object made of two arrays X and Y with npoints each.
Definition TGraph.h:41
Manages Markers.
Definition TMarker.h:22
const char * GetName() const override
Returns name of object.
Definition TNamed.h:47
virtual const char * GetName() const
Returns name of object.
Definition TObject.cxx:439
virtual TObject * Clone(const char *newname="") const
Make a clone of an object using the Streamer facility.
Definition TObject.cxx:223
void Start(Bool_t reset=kTRUE)
Start the stopwatch.
void Stop()
Stop the stopwatch.
void Print(Option_t *option="") const override
Print the real and cpu time passed between the start and stop events.
Basic string class.
Definition TString.h:139
int CompareTo(const char *cs, ECaseCompare cmp=kExact) const
Compare a string to char *cs2.
Definition TString.cxx:450
const Int_t n
Definition legend1.C:16
OwningPtr< T > owningPtr(std::unique_ptr< T > &&ptr)
Internal helper to turn a std::unique_ptr<T> into an OwningPtr.
Definition Config.h:50
T * OwningPtr
An alias for raw pointers for indicating that the return type of a RooFit function is an owning point...
Definition Config.h:43
Definition graph.py:1
Config argument to RooMinimizer ctor.
std::string minimizerType
TMarker m
Definition textangle.C:8