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
21Wrapper 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 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 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{
160 RooSentinel::activate();
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 class 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 _theFitter->FitFCN(*_fcn->getMultiGenFcn());
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 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 }
323 _fcn->Synchronize(_theFitter->Config().ParamsSettings());
324
326 _theFitter->Config().SetMinimizer(_cfg.minimizerType.c_str(), alg);
327
328 profileStart();
329 {
330 auto ctx = makeEvalErrorContext();
331
332 bool ret = fitFcn();
333 _status = ((ret) ? _theFitter->Result().Status() : -1);
334 }
335 profileStop();
336 _fcn->BackProp(_theFitter->Result());
337
338 saveStatus("MINIMIZE", _status);
339
340 return _status;
341}
342
343////////////////////////////////////////////////////////////////////////////////
344/// Execute MIGRAD. Changes in parameter values
345/// and calculated errors are automatically
346/// propagated back the RooRealVars representing
347/// the floating parameters in the MINUIT operation.
348
350{
351 _fcn->Synchronize(_theFitter->Config().ParamsSettings());
352 profileStart();
353 {
354 auto ctx = makeEvalErrorContext();
355
356 _theFitter->Config().SetMinimizer(_cfg.minimizerType.c_str(), "migrad");
357 bool ret = fitFcn();
358 _status = ((ret) ? _theFitter->Result().Status() : -1);
359 }
360 profileStop();
361 _fcn->BackProp(_theFitter->Result());
362
363 saveStatus("MIGRAD", _status);
364
365 return _status;
366}
367
368////////////////////////////////////////////////////////////////////////////////
369/// Execute HESSE. Changes in parameter values
370/// and calculated errors are automatically
371/// propagated back the RooRealVars representing
372/// the floating parameters in the MINUIT operation.
373
375{
376 if (_theFitter->GetMinimizer() == nullptr) {
377 coutW(Minimization) << "RooMinimizer::hesse: Error, run Migrad before Hesse!" << endl;
378 _status = -1;
379 } else {
380
381 _fcn->Synchronize(_theFitter->Config().ParamsSettings());
382 profileStart();
383 {
384 auto ctx = makeEvalErrorContext();
385
386 _theFitter->Config().SetMinimizer(_cfg.minimizerType.c_str());
387 bool ret = _theFitter->CalculateHessErrors();
388 _status = ((ret) ? _theFitter->Result().Status() : -1);
389 }
390 profileStop();
391 _fcn->BackProp(_theFitter->Result());
392
393 saveStatus("HESSE", _status);
394 }
395
396 return _status;
397}
398
399////////////////////////////////////////////////////////////////////////////////
400/// Execute MINOS. Changes in parameter values
401/// and calculated errors are automatically
402/// propagated back the RooRealVars representing
403/// the floating parameters in the MINUIT operation.
404
406{
407 if (_theFitter->GetMinimizer() == nullptr) {
408 coutW(Minimization) << "RooMinimizer::minos: Error, run Migrad before Minos!" << endl;
409 _status = -1;
410 } else {
411
412 _fcn->Synchronize(_theFitter->Config().ParamsSettings());
413 profileStart();
414 {
415 auto ctx = makeEvalErrorContext();
416
417 _theFitter->Config().SetMinimizer(_cfg.minimizerType.c_str());
418 bool ret = _theFitter->CalculateMinosErrors();
419 _status = ((ret) ? _theFitter->Result().Status() : -1);
420 }
421
422 profileStop();
423 _fcn->BackProp(_theFitter->Result());
424
425 saveStatus("MINOS", _status);
426 }
427
428 return _status;
429}
430
431////////////////////////////////////////////////////////////////////////////////
432/// Execute MINOS for given list of parameters. Changes in parameter values
433/// and calculated errors are automatically
434/// propagated back the RooRealVars representing
435/// the floating parameters in the MINUIT operation.
436
437int RooMinimizer::minos(const RooArgSet &minosParamList)
438{
439 if (_theFitter->GetMinimizer() == nullptr) {
440 coutW(Minimization) << "RooMinimizer::minos: Error, run Migrad before Minos!" << endl;
441 _status = -1;
442 } else if (!minosParamList.empty()) {
443
444 _fcn->Synchronize(_theFitter->Config().ParamsSettings());
445 profileStart();
446 {
447 auto ctx = makeEvalErrorContext();
448
449 // get list of parameters for Minos
450 std::vector<unsigned int> paramInd;
451 for (RooAbsArg *arg : minosParamList) {
452 RooAbsArg *par = _fcn->GetFloatParamList()->find(arg->GetName());
453 if (par && !par->isConstant()) {
454 int index = _fcn->GetFloatParamList()->index(par);
455 paramInd.push_back(index);
456 }
457 }
458
459 if (!paramInd.empty()) {
460 // set the parameter indices
461 _theFitter->Config().SetMinosErrors(paramInd);
462
463 _theFitter->Config().SetMinimizer(_cfg.minimizerType.c_str());
464 bool ret = _theFitter->CalculateMinosErrors();
465 _status = ((ret) ? _theFitter->Result().Status() : -1);
466 // to avoid that following minimization computes automatically the Minos errors
467 _theFitter->Config().SetMinosErrors(false);
468 }
469 }
470 profileStop();
471 _fcn->BackProp(_theFitter->Result());
472
473 saveStatus("MINOS", _status);
474 }
475
476 return _status;
477}
478
479////////////////////////////////////////////////////////////////////////////////
480/// Execute SEEK. Changes in parameter values
481/// and calculated errors are automatically
482/// propagated back the RooRealVars representing
483/// the floating parameters in the MINUIT operation.
484
486{
487 _fcn->Synchronize(_theFitter->Config().ParamsSettings());
488 profileStart();
489 {
490 auto ctx = makeEvalErrorContext();
491
492 _theFitter->Config().SetMinimizer(_cfg.minimizerType.c_str(), "seek");
493 bool ret = fitFcn();
494 _status = ((ret) ? _theFitter->Result().Status() : -1);
495 }
496 profileStop();
497 _fcn->BackProp(_theFitter->Result());
498
499 saveStatus("SEEK", _status);
500
501 return _status;
502}
503
504////////////////////////////////////////////////////////////////////////////////
505/// Execute SIMPLEX. Changes in parameter values
506/// and calculated errors are automatically
507/// propagated back the RooRealVars representing
508/// the floating parameters in the MINUIT operation.
509
511{
512 _fcn->Synchronize(_theFitter->Config().ParamsSettings());
513 profileStart();
514 {
515 auto ctx = makeEvalErrorContext();
516
517 _theFitter->Config().SetMinimizer(_cfg.minimizerType.c_str(), "simplex");
518 bool ret = fitFcn();
519 _status = ((ret) ? _theFitter->Result().Status() : -1);
520 }
521 profileStop();
522 _fcn->BackProp(_theFitter->Result());
523
524 saveStatus("SEEK", _status);
525
526 return _status;
527}
528
529////////////////////////////////////////////////////////////////////////////////
530/// Execute IMPROVE. Changes in parameter values
531/// and calculated errors are automatically
532/// propagated back the RooRealVars representing
533/// the floating parameters in the MINUIT operation.
534
536{
537 _fcn->Synchronize(_theFitter->Config().ParamsSettings());
538 profileStart();
539 {
540 auto ctx = makeEvalErrorContext();
541
542 _theFitter->Config().SetMinimizer(_cfg.minimizerType.c_str(), "migradimproved");
543 bool ret = fitFcn();
544 _status = ((ret) ? _theFitter->Result().Status() : -1);
545 }
546 profileStop();
547 _fcn->BackProp(_theFitter->Result());
548
549 saveStatus("IMPROVE", _status);
550
551 return _status;
552}
553
554////////////////////////////////////////////////////////////////////////////////
555/// Change the MINUIT internal printing level
556
558{
559 _theFitter->Config().MinimizerOptions().SetPrintLevel(newLevel + 1);
560}
561
562////////////////////////////////////////////////////////////////////////////////
563/// Get the MINUIT internal printing level
564
566{
567 return _theFitter->Config().MinimizerOptions().PrintLevel() + 1;
568}
569
570////////////////////////////////////////////////////////////////////////////////
571/// If flag is true, perform constant term optimization on
572/// function being minimized.
573
575{
576 _fcn->setOptimizeConst(flag);
577}
578
579////////////////////////////////////////////////////////////////////////////////
580/// Save and return a RooFitResult snapshot of current minimizer status.
581/// This snapshot contains the values of all constant parameters,
582/// the value of all floating parameters at RooMinimizer construction and
583/// after the last MINUIT operation, the MINUIT status, variance quality,
584/// EDM setting, number of calls with evaluation problems, the minimized
585/// function value and the full correlation matrix.
586
587RooFit::OwningPtr<RooFitResult> RooMinimizer::save(const char *userName, const char *userTitle)
588{
589 if (_theFitter->GetMinimizer() == nullptr) {
590 coutW(Minimization) << "RooMinimizer::save: Error, run minimization before!" << endl;
591 return nullptr;
592 }
593
595 TString title;
596 name = userName ? userName : Form("%s", _fcn->getFunctionName().c_str());
597 title = userTitle ? userTitle : Form("%s", _fcn->getFunctionTitle().c_str());
598 auto fitRes = std::make_unique<RooFitResult>(name, title);
599
600 // Move eventual fixed parameters in floatList to constList
601 RooArgList saveConstList(*(_fcn->GetConstParamList()));
602 RooArgList saveFloatInitList(*(_fcn->GetInitFloatParamList()));
603 RooArgList saveFloatFinalList(*(_fcn->GetFloatParamList()));
604 for (std::size_t i = 0; i < _fcn->GetFloatParamList()->size(); i++) {
605 RooAbsArg *par = _fcn->GetFloatParamList()->at(i);
606 if (par->isConstant()) {
607 saveFloatInitList.remove(*saveFloatInitList.find(par->GetName()), true);
608 saveFloatFinalList.remove(*par);
609 saveConstList.add(*par);
610 }
611 }
612 saveConstList.sort();
613
614 fitRes->setConstParList(saveConstList);
615 fitRes->setInitParList(saveFloatInitList);
616
617 double removeOffset = 0.;
618 fitRes->setNumInvalidNLL(_fcn->GetNumInvalidNLL());
619 removeOffset = -_fcn->getOffset();
620
621 fitRes->setStatus(_status);
622 fitRes->setCovQual(_theFitter->GetMinimizer()->CovMatrixStatus());
623 fitRes->setMinNLL(_theFitter->Result().MinFcnValue() + removeOffset);
624 fitRes->setEDM(_theFitter->Result().Edm());
625 fitRes->setFinalParList(saveFloatFinalList);
626 if (!_extV) {
627 std::vector<double> globalCC;
628 TMatrixDSym corrs(_theFitter->Result().Parameters().size());
629 TMatrixDSym covs(_theFitter->Result().Parameters().size());
630 for (std::size_t ic = 0; ic < _theFitter->Result().Parameters().size(); ic++) {
631 globalCC.push_back(_theFitter->Result().GlobalCC(ic));
632 for (std::size_t ii = 0; ii < _theFitter->Result().Parameters().size(); ii++) {
633 corrs(ic, ii) = _theFitter->Result().Correlation(ic, ii);
634 covs(ic, ii) = _theFitter->Result().CovMatrix(ic, ii);
635 }
636 }
637 fitRes->fillCorrMatrix(globalCC, corrs, covs);
638 } else {
639 fitRes->setCovarianceMatrix(*_extV);
640 }
641
642 fitRes->setStatusHistory(_statusHistory);
643
644 return RooFit::Detail::owningPtr(std::move(fitRes));
645}
646
647////////////////////////////////////////////////////////////////////////////////
648/// Create and draw a TH2 with the error contours in the parameters `var1` and `var2`.
649/// \param[in] var1 The first parameter (x axis).
650/// \param[in] var2 The second parameter (y axis).
651/// \param[in] n1 First contour.
652/// \param[in] n2 Optional contour. 0 means don't draw.
653/// \param[in] n3 Optional contour. 0 means don't draw.
654/// \param[in] n4 Optional contour. 0 means don't draw.
655/// \param[in] n5 Optional contour. 0 means don't draw.
656/// \param[in] n6 Optional contour. 0 means don't draw.
657/// \param[in] npoints Number of points for evaluating the contour.
658///
659/// Up to six contours can be drawn using the arguments `n1` to `n6` to request the desired
660/// coverage in units of \f$ \sigma = n^2 \cdot \mathrm{ErrorDef} \f$.
661/// See ROOT::Math::Minimizer::ErrorDef().
662
663RooPlot *RooMinimizer::contour(RooRealVar &var1, RooRealVar &var2, double n1, double n2, double n3, double n4,
664 double n5, double n6, unsigned int npoints)
665{
666 RooArgList *params = _fcn->GetFloatParamList();
667 RooArgList paramSave;
668 params->snapshot(paramSave);
669
670 // Verify that both variables are floating parameters of PDF
671 int index1 = _fcn->GetFloatParamList()->index(&var1);
672 if (index1 < 0) {
673 coutE(Minimization) << "RooMinimizer::contour(" << GetName() << ") ERROR: " << var1.GetName()
674 << " is not a floating parameter of " << _fcn->getFunctionName() << endl;
675 return nullptr;
676 }
677
678 int index2 = _fcn->GetFloatParamList()->index(&var2);
679 if (index2 < 0) {
680 coutE(Minimization) << "RooMinimizer::contour(" << GetName() << ") ERROR: " << var2.GetName()
681 << " is not a floating parameter of PDF " << _fcn->getFunctionName() << endl;
682 return nullptr;
683 }
684
685 // create and draw a frame
686 RooPlot *frame = new RooPlot(var1, var2);
687
688 // draw a point at the current parameter values
689 TMarker *point = new TMarker(var1.getVal(), var2.getVal(), 8);
690 frame->addObject(point);
691
692 // check first if a inimizer is available. If not means
693 // the minimization is not done , so do it
694 if (_theFitter->GetMinimizer() == nullptr) {
695 coutW(Minimization) << "RooMinimizer::contour: Error, run Migrad before contours!" << endl;
696 return frame;
697 }
698
699 // remember our original value of ERRDEF
700 double errdef = _theFitter->GetMinimizer()->ErrorDef();
701
702 double n[6];
703 n[0] = n1;
704 n[1] = n2;
705 n[2] = n3;
706 n[3] = n4;
707 n[4] = n5;
708 n[5] = n6;
709
710 for (int ic = 0; ic < 6; ic++) {
711 if (n[ic] > 0) {
712
713 // set the value corresponding to an n1-sigma contour
714 _theFitter->GetMinimizer()->SetErrorDef(n[ic] * n[ic] * errdef);
715
716 // calculate and draw the contour
717 std::vector<double> xcoor(npoints + 1);
718 std::vector<double> ycoor(npoints + 1);
719 bool ret = _theFitter->GetMinimizer()->Contour(index1, index2, npoints, xcoor.data(), ycoor.data());
720
721 if (!ret) {
722 coutE(Minimization) << "RooMinimizer::contour(" << GetName()
723 << ") ERROR: MINUIT did not return a contour graph for n=" << n[ic] << endl;
724 } else {
725 xcoor[npoints] = xcoor[0];
726 ycoor[npoints] = ycoor[0];
727 TGraph *graph = new TGraph(npoints + 1, xcoor.data(), ycoor.data());
728
729 graph->SetName(Form("contour_%s_n%f", _fcn->getFunctionName().c_str(), n[ic]));
730 graph->SetLineStyle(ic + 1);
731 graph->SetLineWidth(2);
732 graph->SetLineColor(kBlue);
733 frame->addObject(graph, "L");
734 }
735 }
736 }
737
738 // restore the original ERRDEF
739 _theFitter->GetMinimizer()->SetErrorDef(errdef);
740
741 // restore parameter values
742 params->assign(paramSave);
743
744 return frame;
745}
746
747////////////////////////////////////////////////////////////////////////////////
748/// Add parameters in metadata field to process timer
749
751{
752#ifdef ROOFIT_MULTIPROCESS
753 // parameter indices for use in timing heat matrix
754 std::vector<std::string> parameter_names;
755 for (auto &&parameter : *_fcn->GetFloatParamList()) {
756 parameter_names.push_back(parameter->GetName());
757 if (_cfg.verbose) {
758 coutI(Minimization) << "parameter name: " << parameter_names.back() << std::endl;
759 }
760 }
762#else
763 coutI(Minimization) << "Not adding parameters to processtimer because multiprocessing "
764 << "is not enabled." << std::endl;
765#endif
766}
767
768////////////////////////////////////////////////////////////////////////////////
769/// Start profiling timer
770
772{
773 if (_cfg.profile) {
774 _timer.Start();
775 _cumulTimer.Start(_profileStart ? false : true);
776 _profileStart = true;
777 }
778}
779
780////////////////////////////////////////////////////////////////////////////////
781/// Stop profiling timer and report results of last session
782
784{
785 if (_cfg.profile) {
786 _timer.Stop();
788 coutI(Minimization) << "Command timer: ";
789 _timer.Print();
790 coutI(Minimization) << "Session timer: ";
792 }
793}
794
796{
797 auto *fitterFcn = fitter()->GetFCN();
798 return fitterFcn ? fitterFcn : _fcn->getMultiGenFcn();
799}
800
801////////////////////////////////////////////////////////////////////////////////
802/// Apply results of given external covariance matrix. i.e. propagate its errors
803/// to all RRV parameter representations and give this matrix instead of the
804/// HESSE matrix at the next save() call
805
807{
808 _extV.reset(static_cast<TMatrixDSym *>(V.Clone()));
809 _fcn->ApplyCovarianceMatrix(*_extV);
810}
811
813{
815}
816
818{
819 // Import the results of the last fit performed, interpreting
820 // the fit parameters as the given varList of parameters.
821
822 if (_theFitter == nullptr || _theFitter->GetMinimizer() == nullptr) {
823 oocoutE(nullptr, InputArguments) << "RooMinimizer::save: Error, run minimization before!" << endl;
824 return nullptr;
825 }
826
827 // Verify length of supplied varList
828 if (!varList.empty() && varList.size() != _theFitter->Result().NTotalParameters()) {
829 oocoutE(nullptr, InputArguments)
830 << "RooMinimizer::lastMinuitFit: ERROR: supplied variable list must be either empty " << endl
831 << " or match the number of variables of the last fit ("
832 << _theFitter->Result().NTotalParameters() << ")" << endl;
833 return nullptr;
834 }
835
836 // Verify that all members of varList are of type RooRealVar
837 for (RooAbsArg *arg : varList) {
838 if (!dynamic_cast<RooRealVar *>(arg)) {
839 oocoutE(nullptr, InputArguments) << "RooMinimizer::lastMinuitFit: ERROR: variable '" << arg->GetName()
840 << "' is not of type RooRealVar" << endl;
841 return nullptr;
842 }
843 }
844
845 auto res = std::make_unique<RooFitResult>("lastMinuitFit", "Last MINUIT fit");
846
847 // Extract names of fit parameters
848 // and construct corresponding RooRealVars
849 RooArgList constPars("constPars");
850 RooArgList floatPars("floatPars");
851
852 unsigned int i;
853 for (i = 0; i < _theFitter->Result().NTotalParameters(); ++i) {
854
855 TString varName(_theFitter->Result().GetParameterName(i));
856 bool isConst(_theFitter->Result().IsParameterFixed(i));
857
858 double xlo = _theFitter->Config().ParSettings(i).LowerLimit();
859 double xhi = _theFitter->Config().ParSettings(i).UpperLimit();
860 double xerr = _theFitter->Result().Error(i);
861 double xval = _theFitter->Result().Value(i);
862
863 std::unique_ptr<RooRealVar> var;
864 if (varList.empty()) {
865
866 if ((xlo < xhi) && !isConst) {
867 var = std::make_unique<RooRealVar>(varName, varName, xval, xlo, xhi);
868 } else {
869 var = std::make_unique<RooRealVar>(varName, varName, xval);
870 }
871 var->setConstant(isConst);
872 } else {
873
874 var.reset(static_cast<RooRealVar *>(varList.at(i)->Clone()));
875 var->setConstant(isConst);
876 var->setVal(xval);
877 if (xlo < xhi) {
878 var->setRange(xlo, xhi);
879 }
880
881 if (varName.CompareTo(var->GetName())) {
882 oocoutI(nullptr, Eval) << "RooMinimizer::lastMinuitFit: fit parameter '" << varName
883 << "' stored in variable '" << var->GetName() << "'" << endl;
884 }
885 }
886
887 if (isConst) {
888 constPars.addOwned(std::move(var));
889 } else {
890 var->setError(xerr);
891 floatPars.addOwned(std::move(var));
892 }
893 }
894
895 res->setConstParList(constPars);
896 res->setInitParList(floatPars);
897 res->setFinalParList(floatPars);
898 res->setMinNLL(_theFitter->Result().MinFcnValue());
899 res->setEDM(_theFitter->Result().Edm());
900 res->setCovQual(_theFitter->GetMinimizer()->CovMatrixStatus());
901 res->setStatus(_theFitter->Result().Status());
902 std::vector<double> globalCC;
903 TMatrixDSym corrs(_theFitter->Result().Parameters().size());
904 TMatrixDSym covs(_theFitter->Result().Parameters().size());
905 for (unsigned int ic = 0; ic < _theFitter->Result().Parameters().size(); ic++) {
906 globalCC.push_back(_theFitter->Result().GlobalCC(ic));
907 for (unsigned int ii = 0; ii < _theFitter->Result().Parameters().size(); ii++) {
908 corrs(ic, ii) = _theFitter->Result().Correlation(ic, ii);
909 covs(ic, ii) = _theFitter->Result().CovMatrix(ic, ii);
910 }
911 }
912 res->fillCorrMatrix(globalCC, corrs, covs);
913
914 return RooFit::Detail::owningPtr(std::move(res));
915}
916
917/// Try to recover from invalid function values. When invalid function values
918/// are encountered, a penalty term is returned to the minimiser to make it
919/// back off. This sets the strength of this penalty. \note A strength of zero
920/// is equivalent to a constant penalty (= the gradient vanishes, ROOT < 6.24).
921/// Positive values lead to a gradient pointing away from the undefined
922/// regions. Use ~10 to force the minimiser away from invalid function values.
924{
925 _cfg.recoverFromNaN = strength;
926}
927
928bool RooMinimizer::setLogFile(const char *logf)
929{
930 _cfg.logf = logf;
931 if (_cfg.logf) {
932 return _fcn->SetLogFile(_cfg.logf);
933 } else {
934 return false;
935 }
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 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:2468
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:454
Documentation for the abstract class IBaseFunctionMultiDim.
Definition IFunction.h:61
static const std::string & DefaultMinimizerType()
Common abstract base class for objects that represent a value and a "shape" in RooFit.
Definition RooAbsArg.h:79
bool isConstant() const
Check if the "Constant" attribute is set.
Definition RooAbsArg.h:363
TObject * Clone(const char *newname=nullptr) const override
Make a clone of an object using the Streamer facility.
Definition RooAbsArg.h:91
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.
Abstract base class for objects that represent a real value and implements functionality common to al...
Definition RooAbsReal.h:59
double getVal(const RooArgSet *normalisationSet=nullptr) const
Evaluate object.
Definition RooAbsReal.h:103
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
Wrapper class around ROOT::Fit:Fitter that provides a seamless interface between the minimizer functi...
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
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 class from MINUIT (RooMinimizer default 500 * #parameter...
void setStrategy(int istrat)
Change MINUIT strategy to istrat.
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.
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:378
Variable that can be changed from the outside.
Definition RooRealVar.h:37
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:451
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