Logo ROOT  
Reference Guide
RooAbsMinimizerFcn.cxx
Go to the documentation of this file.
1/*****************************************************************************
2 * Project: RooFit *
3 * Package: RooFitCore *
4 * @(#)root/roofitcore:$Id$
5 * Authors: *
6 * AL, Alfio Lazzaro, INFN Milan, alfio.lazzaro@mi.infn.it *
7 * PB, Patrick Bos, Netherlands eScience Center, p.bos@esciencecenter.nl *
8 * *
9 * *
10 * Redistribution and use in source and binary forms, *
11 * with or without modification, are permitted according to the terms *
12 * listed in LICENSE (http://roofit.sourceforge.net/license.txt) *
13 *****************************************************************************/
14
15//////////////////////////////////////////////////////////////////////////////
16//
17// RooAbsMinimizerFcn is an interface class to the ROOT::Math function
18// for minimization. It contains only the "logistics" of synchronizing
19// between Minuit and RooFit. Its subclasses implement actual interfacing
20// to Minuit by subclassing IMultiGenFunction or IMultiGradFunction.
21//
22
23#include "RooAbsMinimizerFcn.h"
24
25#include "RooAbsArg.h"
26#include "RooAbsPdf.h"
27#include "RooArgSet.h"
28#include "RooRealVar.h"
29#include "RooAbsRealLValue.h"
30#include "RooMsgService.h"
31#include "RooMinimizer.h"
32#include "RooNaNPacker.h"
33
34#include "TClass.h"
35#include "TMatrixDSym.h"
36
37#include <fstream>
38#include <iomanip>
39
40using namespace std;
41
43 : _context(context), _verbose(verbose)
44{
45 // Examine parameter list
46 _floatParamList.reset((RooArgList *)paramList.selectByAttrib("Constant", kFALSE));
47 if (_floatParamList->getSize() > 1) {
48 _floatParamList->sort();
49 }
50 _floatParamList->setName("floatParamList");
51
52 _constParamList.reset((RooArgList *)paramList.selectByAttrib("Constant", kTRUE));
53 if (_constParamList->getSize() > 1) {
54 _constParamList->sort();
55 }
56 _constParamList->setName("constParamList");
57
58 // Remove all non-RooRealVar parameters from list (MINUIT cannot handle them)
59 for (unsigned int i = 0; i < _floatParamList->size(); ) { // Note: Counting loop, since removing from collection!
60 const RooAbsArg* arg = (*_floatParamList).at(i);
61 if (!arg->IsA()->InheritsFrom(RooAbsRealLValue::Class())) {
62 oocoutW(_context,Minimization) << "RooAbsMinimizerFcn::RooAbsMinimizerFcn: removing parameter "
63 << arg->GetName() << " from list because it is not of type RooRealVar" << endl;
64 _floatParamList->remove(*arg);
65 } else {
66 ++i;
67 }
68 }
69
70 _nDim = _floatParamList->getSize();
71
72 // Save snapshot of initial lists
75}
76
78 : _context(other._context), _maxFCN(other._maxFCN),
79 _funcOffset(other._funcOffset),
80 _recoverFromNaNStrength(other._recoverFromNaNStrength),
81 _numBadNLL(other._numBadNLL),
82 _printEvalErrors(other._printEvalErrors), _evalCounter(other._evalCounter),
83 _nDim(other._nDim), _optConst(other._optConst),
84 _floatParamList(new RooArgList(*other._floatParamList)), _constParamList(new RooArgList(*other._constParamList)),
85 _initFloatParamList((RooArgList *)other._initFloatParamList->snapshot(kFALSE)),
86 _initConstParamList((RooArgList *)other._initConstParamList->snapshot(kFALSE)),
87 _logfile(other._logfile), _doEvalErrorWall(other._doEvalErrorWall), _verbose(other._verbose)
88{}
89
90
91/// Internal function to synchronize TMinimizer with current
92/// information in RooAbsReal function parameters
93Bool_t RooAbsMinimizerFcn::synchronizeParameterSettings(std::vector<ROOT::Fit::ParameterSettings> &parameters, Bool_t optConst, Bool_t verbose)
94{
95 Bool_t constValChange(kFALSE);
96 Bool_t constStatChange(kFALSE);
97
98 Int_t index(0);
99
100 // Handle eventual migrations from constParamList -> floatParamList
101 for (index = 0; index < _constParamList->getSize(); index++) {
102
103 RooRealVar *par = dynamic_cast<RooRealVar *>(_constParamList->at(index));
104 if (!par)
105 continue;
106
107 RooRealVar *oldpar = dynamic_cast<RooRealVar *>(_initConstParamList->at(index));
108 if (!oldpar)
109 continue;
110
111 // Test if constness changed
112 if (!par->isConstant()) {
113
114 // Remove from constList, add to floatList
115 _constParamList->remove(*par);
116 _floatParamList->add(*par);
117 _initFloatParamList->addClone(*oldpar);
118 _initConstParamList->remove(*oldpar);
119 constStatChange = kTRUE;
120 _nDim++;
121
122 if (verbose) {
124 << "RooAbsMinimizerFcn::synchronize: parameter " << par->GetName() << " is now floating." << endl;
125 }
126 }
127
128 // Test if value changed
129 if (par->getVal() != oldpar->getVal()) {
130 constValChange = kTRUE;
131 if (verbose) {
133 << "RooAbsMinimizerFcn::synchronize: value of constant parameter " << par->GetName() << " changed from "
134 << oldpar->getVal() << " to " << par->getVal() << endl;
135 }
136 }
137 }
138
139 // Update reference list
141
142 // Synchronize MINUIT with function state
143 // Handle floatParamList
144 for (index = 0; index < _floatParamList->getSize(); index++) {
145 RooRealVar *par = dynamic_cast<RooRealVar *>(_floatParamList->at(index));
146
147 if (!par)
148 continue;
149
150 Double_t pstep(0);
151 Double_t pmin(0);
152 Double_t pmax(0);
153
154 if (!par->isConstant()) {
155
156 // Verify that floating parameter is indeed of type RooRealVar
157 if (!par->IsA()->InheritsFrom(RooRealVar::Class())) {
158 oocoutW(_context, Minimization) << "RooAbsMinimizerFcn::fit: Error, non-constant parameter "
159 << par->GetName() << " is not of type RooRealVar, skipping" << endl;
160 _floatParamList->remove(*par);
161 index--;
162 _nDim--;
163 continue;
164 }
165 // make sure the parameter are in dirty state to enable
166 // a real NLL computation when the minimizer calls the function the first time
167 // (see issue #7659)
168 par->setValueDirty();
169
170 // Set the limits, if not infinite
171 if (par->hasMin())
172 pmin = par->getMin();
173 if (par->hasMax())
174 pmax = par->getMax();
175
176 // Calculate step size
177 pstep = par->getError();
178 if (pstep <= 0) {
179 // Floating parameter without error estitimate
180 if (par->hasMin() && par->hasMax()) {
181 pstep = 0.1 * (pmax - pmin);
182
183 // Trim default choice of error if within 2 sigma of limit
184 if (pmax - par->getVal() < 2 * pstep) {
185 pstep = (pmax - par->getVal()) / 2;
186 } else if (par->getVal() - pmin < 2 * pstep) {
187 pstep = (par->getVal() - pmin) / 2;
188 }
189
190 // If trimming results in zero error, restore default
191 if (pstep == 0) {
192 pstep = 0.1 * (pmax - pmin);
193 }
194
195 } else {
196 pstep = 1;
197 }
198 if (verbose) {
200 << "RooAbsMinimizerFcn::synchronize: WARNING: no initial error estimate available for "
201 << par->GetName() << ": using " << pstep << endl;
202 }
203 }
204 } else {
205 pmin = par->getVal();
206 pmax = par->getVal();
207 }
208
209 // new parameter
210 if (index >= Int_t(parameters.size())) {
211
212 if (par->hasMin() && par->hasMax()) {
213 parameters.emplace_back(par->GetName(), par->getVal(), pstep, pmin, pmax);
214 } else {
215 parameters.emplace_back(par->GetName(), par->getVal(), pstep);
216 if (par->hasMin())
217 parameters.back().SetLowerLimit(pmin);
218 else if (par->hasMax())
219 parameters.back().SetUpperLimit(pmax);
220 }
221
222 continue;
223 }
224
225 Bool_t oldFixed = parameters[index].IsFixed();
226 Double_t oldVar = parameters[index].Value();
227 Double_t oldVerr = parameters[index].StepSize();
228 Double_t oldVlo = parameters[index].LowerLimit();
229 Double_t oldVhi = parameters[index].UpperLimit();
230
231 if (par->isConstant() && !oldFixed) {
232
233 // Parameter changes floating -> constant : update only value if necessary
234 if (oldVar != par->getVal()) {
235 parameters[index].SetValue(par->getVal());
236 if (verbose) {
238 << "RooAbsMinimizerFcn::synchronize: value of parameter " << par->GetName() << " changed from "
239 << oldVar << " to " << par->getVal() << endl;
240 }
241 }
242 parameters[index].Fix();
243 constStatChange = kTRUE;
244 if (verbose) {
246 << "RooAbsMinimizerFcn::synchronize: parameter " << par->GetName() << " is now fixed." << endl;
247 }
248
249 } else if (par->isConstant() && oldFixed) {
250
251 // Parameter changes constant -> constant : update only value if necessary
252 if (oldVar != par->getVal()) {
253 parameters[index].SetValue(par->getVal());
254 constValChange = kTRUE;
255
256 if (verbose) {
258 << "RooAbsMinimizerFcn::synchronize: value of fixed parameter " << par->GetName() << " changed from "
259 << oldVar << " to " << par->getVal() << endl;
260 }
261 }
262
263 } else {
264 // Parameter changes constant -> floating
265 if (!par->isConstant() && oldFixed) {
266 parameters[index].Release();
267 constStatChange = kTRUE;
268
269 if (verbose) {
271 << "RooAbsMinimizerFcn::synchronize: parameter " << par->GetName() << " is now floating." << endl;
272 }
273 }
274
275 // Parameter changes constant -> floating : update all if necessary
276 if (oldVar != par->getVal() || oldVlo != pmin || oldVhi != pmax || oldVerr != pstep) {
277 parameters[index].SetValue(par->getVal());
278 parameters[index].SetStepSize(pstep);
279 if (par->hasMin() && par->hasMax())
280 parameters[index].SetLimits(pmin, pmax);
281 else if (par->hasMin())
282 parameters[index].SetLowerLimit(pmin);
283 else if (par->hasMax())
284 parameters[index].SetUpperLimit(pmax);
285 }
286
287 // Inform user about changes in verbose mode
288 if (verbose) {
289 // if ierr<0, par was moved from the const list and a message was already printed
290
291 if (oldVar != par->getVal()) {
293 << "RooAbsMinimizerFcn::synchronize: value of parameter " << par->GetName() << " changed from "
294 << oldVar << " to " << par->getVal() << endl;
295 }
296 if (oldVlo != pmin || oldVhi != pmax) {
298 << "RooAbsMinimizerFcn::synchronize: limits of parameter " << par->GetName() << " changed from ["
299 << oldVlo << "," << oldVhi << "] to [" << pmin << "," << pmax << "]" << endl;
300 }
301
302 // If oldVerr=0, then parameter was previously fixed
303 if (oldVerr != pstep && oldVerr != 0) {
305 << "RooAbsMinimizerFcn::synchronize: error/step size of parameter " << par->GetName()
306 << " changed from " << oldVerr << " to " << pstep << endl;
307 }
308 }
309 }
310 }
311
312 if (optConst) {
313 optimizeConstantTerms(constStatChange, constValChange);
314 }
315
316 return 0;
317}
318
319Bool_t
320RooAbsMinimizerFcn::Synchronize(std::vector<ROOT::Fit::ParameterSettings> &parameters, Bool_t optConst, Bool_t verbose) {
321 return synchronizeParameterSettings(parameters, optConst, verbose);
322}
323
324/// Modify PDF parameter error by ordinal index (needed by MINUIT)
326{
327 static_cast<RooRealVar*>(_floatParamList->at(index))->setError(value);
328}
329
330/// Modify PDF parameter error by ordinal index (needed by MINUIT)
332{
333 static_cast<RooRealVar*>(_floatParamList->at(index))->removeAsymError();
334}
335
336/// Modify PDF parameter error by ordinal index (needed by MINUIT)
338{
339 static_cast<RooRealVar*>(_floatParamList->at(index))->setAsymError(loVal, hiVal);
340}
341
342/// Transfer MINUIT fit results back into RooFit objects.
344{
345 for (unsigned int index = 0; index < _nDim; index++) {
346 Double_t value = results.Value(index);
347 SetPdfParamVal(index, value);
348
349 // Set the parabolic error
350 Double_t err = results.Error(index);
351 SetPdfParamErr(index, err);
352
353 Double_t eminus = results.LowerError(index);
354 Double_t eplus = results.UpperError(index);
355
356 if (eplus > 0 || eminus < 0) {
357 // Store the asymmetric error, if it is available
358 SetPdfParamErr(index, eminus, eplus);
359 } else {
360 // Clear the asymmetric error
362 }
363 }
364}
365
366/// Change the file name for logging of a RooMinimizer of all MINUIT steppings
367/// through the parameter space. If inLogfile is null, the current log file
368/// is closed and logging is stopped.
370{
371 if (_logfile) {
372 oocoutI(_context, Minimization) << "RooAbsMinimizerFcn::setLogFile: closing previous log file" << endl;
373 _logfile->close();
374 delete _logfile;
375 _logfile = 0;
376 }
377 _logfile = new ofstream(inLogfile);
378 if (!_logfile->good()) {
379 oocoutI(_context, Minimization) << "RooAbsMinimizerFcn::setLogFile: cannot open file " << inLogfile << endl;
380 _logfile->close();
381 delete _logfile;
382 _logfile = 0;
383 }
384
385 return kFALSE;
386}
387
388/// Apply results of given external covariance matrix. i.e. propagate its errors
389/// to all RRV parameter representations and give this matrix instead of the
390/// HESSE matrix at the next save() call
392{
393 for (unsigned int i = 0; i < _nDim; i++) {
394 // Skip fixed parameters
395 if (_floatParamList->at(i)->isConstant()) {
396 continue;
397 }
398 SetPdfParamErr(i, sqrt(V(i, i)));
399 }
400}
401
402/// Set value of parameter i.
403Bool_t RooAbsMinimizerFcn::SetPdfParamVal(int index, double value) const
404{
405 auto par = static_cast<RooRealVar*>(&(*_floatParamList)[index]);
406
407 if (par->getVal()!=value) {
408 if (_verbose) cout << par->GetName() << "=" << value << ", " ;
409
410 par->setVal(value);
411 return kTRUE;
412 }
413
414 return kFALSE;
415}
416
417
418/// Print information about why evaluation failed.
419/// Using _printEvalErrors, the number of errors printed can be steered.
420/// Negative values disable printing.
422 if (_printEvalErrors < 0)
423 return;
424
425 std::ostringstream msg;
426 if (_doEvalErrorWall) {
427 msg << "RooAbsMinimizerFcn: Minimized function has error status." << endl
428 << "Returning maximum FCN so far (" << _maxFCN
429 << ") to force MIGRAD to back out of this region. Error log follows.\n";
430 } else {
431 msg << "RooAbsMinimizerFcn: Minimized function has error status but is ignored.\n";
432 }
433
434 msg << "Parameter values: " ;
435 for (const auto par : *_floatParamList) {
436 auto var = static_cast<const RooRealVar*>(par);
437 msg << "\t" << var->GetName() << "=" << var->getVal() ;
438 }
439 msg << std::endl;
440
442 ooccoutW(_context,Minimization) << msg.str() << endl;
443}
444
445
446////////////////////////////////////////////////////////////////////////////////
447/// Logistics
448
450{
451 return _floatParamList.get();
452}
454{
455 return _constParamList.get();
456}
458{
459 return _initFloatParamList.get();
460}
462{
463 return _initConstParamList.get();
464}
465
467{
468 _doEvalErrorWall = flag;
469}
471{
472 _printEvalErrors = numEvalErrors;
473}
474
476{
477 return _maxFCN;
478}
480{
481 return _numBadNLL;
482}
483
485{
486 return _evalCounter;
487}
489{
490 _evalCounter = 0;
491}
492
494{
495 _verbose = flag;
496}
497
499{
501
502 if (_optConst && !flag) {
503 if (_context->getPrintLevel() > -1)
504 oocoutI(_context, Minimization) << "RooAbsMinimizerFcn::setOptimizeConst: deactivating const optimization" << endl;
506 _optConst = flag;
507 } else if (!_optConst && flag) {
508 if (_context->getPrintLevel() > -1)
509 oocoutI(_context, Minimization) << "RooAbsMinimizerFcn::setOptimizeConst: activating const optimization" << endl;
511 _optConst = flag;
512 } else if (_optConst && flag) {
513 if (_context->getPrintLevel() > -1)
514 oocoutI(_context, Minimization) << "RooAbsMinimizerFcn::setOptimizeConst: const optimization already active" << endl;
515 } else {
516 if (_context->getPrintLevel() > -1)
517 oocoutI(_context, Minimization) << "RooAbsMinimizerFcn::setOptimizeConst: const optimization wasn't active" << endl;
518 }
519
521}
522
523void RooAbsMinimizerFcn::optimizeConstantTerms(bool constStatChange, bool constValChange) {
524 if (constStatChange) {
525
527
528 oocoutI(_context,Minimization) << "RooAbsMinimizerFcn::optimizeConstantTerms: set of constant parameters changed, rerunning const optimizer" << endl ;
530 } else if (constValChange) {
531 oocoutI(_context,Minimization) << "RooAbsMinimizerFcn::optimizeConstantTerms: constant parameter values changed, rerunning const optimizer" << endl ;
533 }
534
536}
537
539{
540 return _optConst;
541}
542
543std::vector<double> RooAbsMinimizerFcn::getParameterValues() const
544{
545 // TODO: make a cache for this somewhere so it doesn't have to be recreated on each call
546 std::vector<double> values;
547 values.reserve(_nDim);
548
549 for (std::size_t index = 0; index < _nDim; ++index) {
550 RooRealVar *par = (RooRealVar *)_floatParamList->at(index);
551 values.push_back(par->getVal());
552 }
553
554 return values;
555}
void Class()
Definition: Class.C:29
#define oocoutW(o, a)
Definition: RooMsgService.h:47
#define oocoutI(o, a)
Definition: RooMsgService.h:45
#define ooccoutW(o, a)
Definition: RooMsgService.h:55
int Int_t
Definition: RtypesCore.h:45
const Bool_t kFALSE
Definition: RtypesCore.h:101
bool Bool_t
Definition: RtypesCore.h:63
const Bool_t kTRUE
Definition: RtypesCore.h:100
class containg the result of the fit and all the related information (fitted parameter values,...
Definition: FitResult.h:47
double UpperError(unsigned int i) const
upper Minos error. If Minos has not run for parameter i return the parabolic error
Definition: FitResult.cxx:383
double Error(unsigned int i) const
parameter error by index
Definition: FitResult.h:186
double Value(unsigned int i) const
parameter value by index
Definition: FitResult.h:179
double LowerError(unsigned int i) const
lower Minos error. If Minos has not run for parameter i return the parabolic error
Definition: FitResult.cxx:376
RooAbsArg is the common abstract base class for objects that represent a value and a "shape" in RooFi...
Definition: RooAbsArg.h:75
void setValueDirty()
Mark the element dirty. This forces a re-evaluation when a value is requested.
Definition: RooAbsArg.h:511
@ DeActivate
Definition: RooAbsArg.h:409
@ ValueChange
Definition: RooAbsArg.h:409
@ ConfigChange
Definition: RooAbsArg.h:409
Bool_t isConstant() const
Check if the "Constant" attribute is set.
Definition: RooAbsArg.h:383
RooAbsCollection * selectByAttrib(const char *name, Bool_t value) const
Create a subset of the current collection, consisting only of those elements with the specified attri...
RooArgList * GetInitFloatParamList()
Bool_t SetPdfParamVal(int index, double value) const
Set value of parameter i.
void setOptimizeConst(Int_t flag)
virtual Bool_t Synchronize(std::vector< ROOT::Fit::ParameterSettings > &parameters, Bool_t optConst, Bool_t verbose)
Like synchronizeParameterSettings, Synchronize informs Minuit through its parameter_settings vector o...
RooMinimizer * _context
std::unique_ptr< RooArgList > _floatParamList
void SetPrintEvalErrors(Int_t numEvalErrors)
Int_t GetNumInvalidNLL() const
RooArgList * GetConstParamList()
virtual void setOptimizeConstOnFunction(RooAbsArg::ConstOpCode opcode, Bool_t doAlsoTrackingOpt)=0
This function must be overridden in the derived class to pass on constant term optimization configura...
void BackProp(const ROOT::Fit::FitResult &results)
Put Minuit results back into RooFit objects.
void optimizeConstantTerms(bool constStatChange, bool constValChange)
void printEvalErrors() const
Print information about why evaluation failed.
void ClearPdfParamAsymErr(Int_t index)
Modify PDF parameter error by ordinal index (needed by MINUIT)
std::vector< double > getParameterValues() const
std::ofstream * _logfile
void ApplyCovarianceMatrix(TMatrixDSym &V)
Set different external covariance matrix.
std::unique_ptr< RooArgList > _initConstParamList
Bool_t synchronizeParameterSettings(std::vector< ROOT::Fit::ParameterSettings > &parameters, Bool_t optConst, Bool_t verbose)
Informs Minuit through its parameter_settings vector of RooFit parameter properties.
RooArgList * GetFloatParamList()
Logistics.
void SetEvalErrorWall(Bool_t flag)
void SetVerbose(Bool_t flag=kTRUE)
std::unique_ptr< RooArgList > _constParamList
RooAbsMinimizerFcn(RooArgList paramList, RooMinimizer *context, bool verbose=false)
RooArgList * GetInitConstParamList()
std::unique_ptr< RooArgList > _initFloatParamList
Bool_t SetLogFile(const char *inLogfile)
Change the file name for logging of a RooMinimizer of all MINUIT steppings through the parameter spac...
void SetPdfParamErr(Int_t index, Double_t value)
Modify PDF parameter error by ordinal index (needed by MINUIT)
virtual Double_t getMax(const char *name=0) const
Get maximum of currently defined range.
Bool_t hasMax(const char *name=0) const
Check if variable has an upper bound.
Bool_t hasMin(const char *name=0) const
Check if variable has a lower bound.
virtual Double_t getMin(const char *name=0) const
Get miniminum of currently defined range.
Double_t getVal(const RooArgSet *normalisationSet=nullptr) const
Evaluate object.
Definition: RooAbsReal.h:93
static void setEvalErrorLoggingMode(ErrorLoggingMode m)
Set evaluation error logging mode.
static void printEvalErrors(std::ostream &os=std::cout, Int_t maxPerNode=10000000)
Print all outstanding logged evaluation error on the given ostream.
RooArgList is a container object that can hold multiple RooAbsArg objects.
Definition: RooArgList.h:22
RooMinimizer is a wrapper class around ROOT::Fit:Fitter that provides a seamless interface between th...
Definition: RooMinimizer.h:49
Int_t getPrintLevel() const
RooRealVar represents a variable that can be changed from the outside.
Definition: RooRealVar.h:39
Double_t getError() const
Definition: RooRealVar.h:60
virtual const char * GetName() const
Returns name of object.
Definition: TNamed.h:47
virtual Bool_t InheritsFrom(const char *classname) const
Returns kTRUE if object inherits from class "classname".
Definition: TObject.cxx:445
VecExpr< UnaryOp< Sqrt< T >, VecExpr< A, T, D >, T >, T, D > sqrt(const VecExpr< A, T, D > &rhs)
@ Minimization
Definition: RooGlobalFunc.h:60
static constexpr double eplus