Logo ROOT  
Reference Guide
RooFitResult.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 * *
9 * Copyright (c) 2000-2005, Regents of the University of California *
10 * and Stanford University. All rights reserved. *
11 * *
12 * Redistribution and use in source and binary forms, *
13 * with or without modification, are permitted according to the terms *
14 * listed in LICENSE (http://roofit.sourceforge.net/license.txt) *
15 *****************************************************************************/
16
17//////////////////////////////////////////////////////////////////////////////
18/// \class RooFitResult
19/// RooFitResult is a container class to hold the input and output
20/// of a PDF fit to a dataset. It contains:
21///
22/// * Values of all constant parameters
23/// * Initial and final values of floating parameters with error
24/// * Correlation matrix and global correlation coefficients
25/// * NLL and EDM at mininum
26///
27/// No references to the fitted PDF and dataset are stored
28///
29
30#include <iostream>
31#include <iomanip>
32
33#include "TBuffer.h"
34#include "TMinuit.h"
35#include "TMath.h"
36#include "TMarker.h"
37#include "TLine.h"
38#include "TBox.h"
39#include "TGaxis.h"
40#include "TMatrix.h"
41#include "TVector.h"
42#include "TDirectory.h"
43#include "TClass.h"
44#include "RooFitResult.h"
45#include "RooArgSet.h"
46#include "RooArgList.h"
47#include "RooRealVar.h"
48#include "RooPlot.h"
49#include "RooEllipse.h"
50#include "RooRandom.h"
51#include "RooMsgService.h"
52#include "TH2D.h"
53#include "TMatrixDSym.h"
54#include "RooMultiVarGaussian.h"
55
56
57using namespace std;
58
60
61
62
63////////////////////////////////////////////////////////////////////////////////
64/// Constructor with name and title
65
66RooFitResult::RooFitResult(const char* name, const char* title) :
67 TNamed(name,title), _constPars(0), _initPars(0), _finalPars(0), _globalCorr(0), _randomPars(0), _Lt(0),
68 _CM(0), _VM(0), _GC(0)
69{
70 if (name) appendToDir(this,true) ;
71}
72
73
74////////////////////////////////////////////////////////////////////////////////
75/// Copy constructor
76
78 TNamed(other),
79 RooPrintable(other),
80 RooDirItem(other),
81 _status(other._status),
82 _covQual(other._covQual),
83 _numBadNLL(other._numBadNLL),
84 _minNLL(other._minNLL),
85 _edm(other._edm),
86 _globalCorr(0),
87 _randomPars(0),
88 _Lt(0),
89 _CM(0),
90 _VM(0),
91 _GC(0),
92 _statusHistory(other._statusHistory)
93{
100 if (other._randomPars) {
103 }
104 if (other._Lt) _Lt = new TMatrix(*other._Lt);
105 if (other._VM) _VM = new TMatrixDSym(*other._VM) ;
106 if (other._CM) _CM = new TMatrixDSym(*other._CM) ;
107 if (other._GC) _GC = new TVectorD(*other._GC) ;
108
109 if (GetName())
110 appendToDir(this, true);
111}
112
113
114
115////////////////////////////////////////////////////////////////////////////////
116/// Destructor
117
119{
120 if (_constPars) delete _constPars ;
121 if (_initPars) delete _initPars ;
122 if (_finalPars) delete _finalPars ;
123 if (_globalCorr) delete _globalCorr;
124 if (_randomPars) delete _randomPars;
125 if (_Lt) delete _Lt;
126 if (_CM) delete _CM ;
127 if (_VM) delete _VM ;
128 if (_GC) delete _GC ;
129
132
133 removeFromDir(this) ;
134}
135
136
137////////////////////////////////////////////////////////////////////////////////
138/// Fill the list of constant parameters
139
141{
142 if (_constPars) delete _constPars ;
144 list.snapshot(*_constPars);
145 for(auto* rrv : dynamic_range_cast<RooRealVar*>(*_constPars)) {
146 if (rrv) {
147 rrv->deleteSharedProperties() ;
148 }
149 }
150}
151
152
153
154////////////////////////////////////////////////////////////////////////////////
155/// Fill the list of initial values of the floating parameters
156
158{
159 if (_initPars) delete _initPars ;
160 _initPars = new RooArgList;
161 list.snapshot(*_initPars);
162 for(auto* rrv : dynamic_range_cast<RooRealVar*>(*_initPars)) {
163 if (rrv) {
164 rrv->deleteSharedProperties() ;
165 }
166 }
167}
168
169
170
171////////////////////////////////////////////////////////////////////////////////
172/// Fill the list of final values of the floating parameters
173
175{
176 if (_finalPars) delete _finalPars ;
178 list.snapshot(*_finalPars);
179
180 for(auto* rrv : dynamic_range_cast<RooRealVar*>(*_finalPars)) {
181 if (rrv) {
182 rrv->deleteSharedProperties() ;
183 }
184 }
185}
186
187
188
189////////////////////////////////////////////////////////////////////////////////
190
192{
193 if (icycle>=_statusHistory.size()) {
194 coutE(InputArguments) << "RooFitResult::statusCodeHistory(" << GetName()
195 << " ERROR request for status history slot "
196 << icycle << " exceeds history count of " << _statusHistory.size() << endl ;
197 }
198 return _statusHistory[icycle].second ;
199}
200
201
202
203////////////////////////////////////////////////////////////////////////////////
204
206{
207 if (icycle>=_statusHistory.size()) {
208 coutE(InputArguments) << "RooFitResult::statusLabelHistory(" << GetName()
209 << " ERROR request for status history slot "
210 << icycle << " exceeds history count of " << _statusHistory.size() << endl ;
211 }
212 return _statusHistory[icycle].first.c_str() ;
213}
214
215
216
217////////////////////////////////////////////////////////////////////////////////
218/// Add objects to a 2D plot that represent the fit results for the
219/// two named parameters. The input frame with the objects added is
220/// returned, or zero in case of an error. Which objects are added
221/// are determined by the options string which should be a concatenation
222/// of the following (not case sensitive):
223///
224/// * M - a marker at the best fit result
225/// * E - an error ellipse calculated at 1-sigma using the error matrix at the minimum
226/// * 1 - the 1-sigma error bar for parameter 1
227/// * 2 - the 1-sigma error bar for parameter 2
228/// * B - the bounding box for the error ellipse
229/// * H - a line and horizontal axis for reading off the correlation coefficient
230/// * V - a line and vertical axis for reading off the correlation coefficient
231/// * A - draw axes for reading off the correlation coefficients with the H or V options
232///
233/// You can change the attributes of objects in the returned RooPlot using the
234/// various `RooPlot::getAttXxx(name)` member functions, e.g.
235/// ```
236/// plot->getAttLine("contour")->SetLineStyle(kDashed);
237/// ```
238/// Use `plot->Print()` for a list of all objects and their names (unfortunately most
239/// of the ROOT builtin graphics objects like TLine are unnamed). Drag the left mouse
240/// button along the labels of either axis button to interactively zoom in a plot.
241
242RooPlot *RooFitResult::plotOn(RooPlot *frame, const char *parName1, const char *parName2,
243 const char *options) const
244{
245 // lookup the input parameters by name: we require that they were floated in our fit
246 const RooRealVar *par1= dynamic_cast<const RooRealVar*>(floatParsFinal().find(parName1));
247 if(0 == par1) {
248 coutE(InputArguments) << "RooFitResult::correlationPlot: parameter not floated in fit: " << parName1 << endl;
249 return 0;
250 }
251 const RooRealVar *par2= dynamic_cast<const RooRealVar*>(floatParsFinal().find(parName2));
252 if(0 == par2) {
253 coutE(InputArguments) << "RooFitResult::correlationPlot: parameter not floated in fit: " << parName2 << endl;
254 return 0;
255 }
256
257 // options are not case sensitive
258 TString opt(options);
259 opt.ToUpper();
260
261 // lookup the 2x2 covariance matrix elements for these variables
262 double x1= par1->getVal();
263 double x2= par2->getVal();
264 double s1= par1->getError();
265 double s2= par2->getError();
266 double rho= correlation(parName1, parName2);
267
268 // add a 1-sigma error ellipse, if requested
269 if(opt.Contains("E")) {
270 RooEllipse *contour= new RooEllipse("contour",x1,x2,s1,s2,rho);
271 contour->SetLineWidth(2) ;
272 frame->addPlotable(contour);
273 }
274
275 // add the error bar for parameter 1, if requested
276 if(opt.Contains("1")) {
277 TLine *hline= new TLine(x1-s1,x2,x1+s1,x2);
278 hline->SetLineColor(kRed);
279 frame->addObject(hline);
280 }
281
282 if(opt.Contains("2")) {
283 TLine *vline= new TLine(x1,x2-s2,x1,x2+s2);
284 vline->SetLineColor(kRed);
285 frame->addObject(vline);
286 }
287
288 if(opt.Contains("B")) {
289 TBox *box= new TBox(x1-s1,x2-s2,x1+s1,x2+s2);
290 box->SetLineStyle(kDashed);
291 box->SetLineColor(kRed);
292 box->SetFillStyle(0);
293 frame->addObject(box);
294 }
295
296 if(opt.Contains("H")) {
297 TLine *line= new TLine(x1-rho*s1,x2-s2,x1+rho*s1,x2+s2);
300 line->SetLineWidth(2) ;
301 frame->addObject(line);
302 if(opt.Contains("A")) {
303 TGaxis *axis= new TGaxis(x1-s1,x2-s2,x1+s1,x2-s2,-1.,+1.,502,"-=");
304 axis->SetLineColor(kBlue);
305 frame->addObject(axis);
306 }
307 }
308
309 if(opt.Contains("V")) {
310 TLine *line= new TLine(x1-s1,x2-rho*s2,x1+s1,x2+rho*s2);
313 line->SetLineWidth(2) ;
314 frame->addObject(line);
315 if(opt.Contains("A")) {
316 TGaxis *axis= new TGaxis(x1-s1,x2-s2,x1-s1,x2+s2,-1.,+1.,502,"-=");
317 axis->SetLineColor(kBlue);
318 frame->addObject(axis);
319 }
320 }
321
322 // add a marker at the fitted value, if requested
323 if(opt.Contains("M")) {
324 TMarker *marker= new TMarker(x1,x2,20);
325 marker->SetMarkerColor(kBlack);
326 frame->addObject(marker);
327 }
328
329 return frame;
330}
331
332
333////////////////////////////////////////////////////////////////////////////////
334/// Return a list of floating parameter values that are perturbed from the final
335/// fit values by random amounts sampled from the covariance matrix. The returned
336/// object is overwritten with each call and belongs to the RooFitResult. Uses
337/// the "square root method" to decompose the covariance matrix, which makes inverting
338/// it unnecessary.
339
341{
342 Int_t nPar= _finalPars->getSize();
343 if(0 == _randomPars) { // first-time initialization
344 assert(0 != _finalPars);
345 // create the list of random values to fill
348 // calculate the elements of the upper-triangular matrix L that gives Lt*L = C
349 // where Lt is the transpose of L (the "square-root method")
350 TMatrix L(nPar,nPar);
351 for(Int_t iPar= 0; iPar < nPar; iPar++) {
352 // calculate the diagonal term first
353 L(iPar,iPar)= covariance(iPar,iPar);
354 for(Int_t k= 0; k < iPar; k++) {
355 double tmp= L(k,iPar);
356 L(iPar,iPar)-= tmp*tmp;
357 }
358 L(iPar,iPar)= sqrt(L(iPar,iPar));
359 // then the off-diagonal terms
360 for(Int_t jPar= iPar+1; jPar < nPar; jPar++) {
361 L(iPar,jPar)= covariance(iPar,jPar);
362 for(Int_t k= 0; k < iPar; k++) {
363 L(iPar,jPar)-= L(k,iPar)*L(k,jPar);
364 }
365 L(iPar,jPar)/= L(iPar,iPar);
366 }
367 }
368 // remember Lt
370 }
371 else {
372 // reset to the final fit values
374 }
375
376 // create a vector of unit Gaussian variables
377 TVector g(nPar);
378 for(Int_t k= 0; k < nPar; k++) g(k)= RooRandom::gaussian();
379 // multiply this vector by Lt to introduce the appropriate correlations
380 g*= (*_Lt);
381 // add the mean value offsets and store the results
382 Int_t index(0);
383 for(auto * par : static_range_cast<RooRealVar*>(*_randomPars)) {
384 par->setVal(par->getVal() + g(index++));
385 }
386
387 return *_randomPars;
388}
389
390
391////////////////////////////////////////////////////////////////////////////////
392/// Return the correlation between parameters 'par1' and 'par2'
393
394double RooFitResult::correlation(const char* parname1, const char* parname2) const
395{
396 Int_t idx1 = _finalPars->index(parname1) ;
397 Int_t idx2 = _finalPars->index(parname2) ;
398 if (idx1<0) {
399 coutE(InputArguments) << "RooFitResult::correlation(" << GetName() << ") parameter " << parname1 << " is not a floating fit parameter" << endl ;
400 return 0 ;
401 }
402 if (idx2<0) {
403 coutE(InputArguments) << "RooFitResult::correlation(" << GetName() << ") parameter " << parname2 << " is not a floating fit parameter" << endl ;
404 return 0 ;
405 }
406 return correlation(idx1,idx2) ;
407}
408
409
410
411////////////////////////////////////////////////////////////////////////////////
412/// Return the set of correlation coefficients of parameter 'par' with
413/// all other floating parameters
414
415const RooArgList* RooFitResult::correlation(const char* parname) const
416{
417 if (_globalCorr==0) {
419 }
420
421 RooAbsArg* arg = _initPars->find(parname) ;
422 if (!arg) {
423 coutE(InputArguments) << "RooFitResult::correlation: variable " << parname << " not a floating parameter in fit" << endl ;
424 return 0 ;
425 }
426 return (RooArgList*)_corrMatrix.At(_initPars->index(arg)) ;
427}
428
429
430
431////////////////////////////////////////////////////////////////////////////////
432/// Return the global correlation of the named parameter
433
434double RooFitResult::globalCorr(const char* parname)
435{
436 if (_globalCorr==0) {
438 }
439
440 RooAbsArg* arg = _initPars->find(parname) ;
441 if (!arg) {
442 coutE(InputArguments) << "RooFitResult::globalCorr: variable " << parname << " not a floating parameter in fit" << endl ;
443 return 0 ;
444 }
445
446 if (_globalCorr) {
447 return ((RooAbsReal*)_globalCorr->at(_initPars->index(arg)))->getVal() ;
448 } else {
449 return 1.0 ;
450 }
451}
452
453
454
455////////////////////////////////////////////////////////////////////////////////
456/// Return the list of all global correlations
457
459{
460 if (_globalCorr==0) {
462 }
463
464 return _globalCorr ;
465}
466
467
468
469////////////////////////////////////////////////////////////////////////////////
470/// Return a correlation matrix element addressed with numeric indices.
471
473{
474 return (*_CM)(row,col) ;
475}
476
477
478////////////////////////////////////////////////////////////////////////////////
479/// Return the covariance matrix element addressed with numeric indices.
480
482{
483 return (*_VM)(row,col) ;
484}
485
486
487
488////////////////////////////////////////////////////////////////////////////////
489/// Print fit result to stream 'os'. In Verbose mode, the constant parameters and
490/// the initial and final values of the floating parameters are printed.
491/// Standard mode only the final values of the floating parameters are printed
492
493void RooFitResult::printMultiline(ostream& os, Int_t /*contents*/, bool verbose, TString indent) const
494{
495
496 os << endl
497 << indent << " RooFitResult: minimized FCN value: " << _minNLL << ", estimated distance to minimum: " << _edm << endl
498 << indent << " covariance matrix quality: " ;
499 switch(_covQual) {
500 case -1 : os << "Unknown, matrix was externally provided" ; break ;
501 case 0 : os << "Not calculated at all" ; break ;
502 case 1 : os << "Approximation only, not accurate" ; break ;
503 case 2 : os << "Full matrix, but forced positive-definite" ; break ;
504 case 3 : os << "Full, accurate covariance matrix" ; break ;
505 }
506 os << endl ;
507 os << indent << " Status : " ;
508 for (vector<pair<string,int> >::const_iterator iter = _statusHistory.begin() ; iter != _statusHistory.end() ; ++iter) {
509 os << iter->first << "=" << iter->second << " " ;
510 }
511 os << endl << endl ;;
512
513 Int_t i ;
514 if (verbose) {
515 if (_constPars->getSize()>0) {
516 os << indent << " Constant Parameter Value " << endl
517 << indent << " -------------------- ------------" << endl ;
518
519 for (i=0 ; i<_constPars->getSize() ; i++) {
520 os << indent << " " << setw(20) << ((RooAbsArg*)_constPars->at(i))->GetName()
521 << " " << setw(12) << Form("%12.4e",((RooRealVar*)_constPars->at(i))->getVal())
522 << endl ;
523 }
524
525 os << endl ;
526 }
527
528 // Has any parameter asymmetric errors?
529 bool doAsymErr(false) ;
530 for (i=0 ; i<_finalPars->getSize() ; i++) {
531 if (((RooRealVar*)_finalPars->at(i))->hasAsymError()) {
532 doAsymErr=true ;
533 break ;
534 }
535 }
536
537 if (doAsymErr) {
538 os << indent << " Floating Parameter InitialValue FinalValue (+HiError,-LoError) GblCorr." << endl
539 << indent << " -------------------- ------------ ---------------------------------- --------" << endl ;
540 } else {
541 os << indent << " Floating Parameter InitialValue FinalValue +/- Error GblCorr." << endl
542 << indent << " -------------------- ------------ -------------------------- --------" << endl ;
543 }
544
545 for (i=0 ; i<_finalPars->getSize() ; i++) {
546 os << indent << " " << setw(20) << ((RooAbsArg*)_finalPars->at(i))->GetName() ;
547 os << indent << " " << setw(12) << Form("%12.4e",((RooRealVar*)_initPars->at(i))->getVal())
548 << indent << " " << setw(12) << Form("%12.4e",((RooRealVar*)_finalPars->at(i))->getVal()) ;
549
550 if (((RooRealVar*)_finalPars->at(i))->hasAsymError()) {
551 os << setw(21) << Form(" (+%8.2e,-%8.2e)",((RooRealVar*)_finalPars->at(i))->getAsymErrorHi(),
552 -1*((RooRealVar*)_finalPars->at(i))->getAsymErrorLo()) ;
553 } else {
554 double err = ((RooRealVar*)_finalPars->at(i))->getError() ;
555 os << (doAsymErr?" ":"") << " +/- " << setw(9) << Form("%9.2e",err) ;
556 }
557
558 if (_globalCorr) {
559 os << " " << setw(8) << Form("%8.6f" ,((RooRealVar*)_globalCorr->at(i))->getVal()) ;
560 } else {
561 os << " <none>" ;
562 }
563
564 os << endl ;
565 }
566
567 } else {
568 os << indent << " Floating Parameter FinalValue +/- Error " << endl
569 << indent << " -------------------- --------------------------" << endl ;
570
571 for (i=0 ; i<_finalPars->getSize() ; i++) {
572 double err = ((RooRealVar*)_finalPars->at(i))->getError() ;
573 os << indent << " " << setw(20) << ((RooAbsArg*)_finalPars->at(i))->GetName()
574 << " " << setw(12) << Form("%12.4e",((RooRealVar*)_finalPars->at(i))->getVal())
575 << " +/- " << setw(9) << Form("%9.2e",err)
576 << endl ;
577 }
578 }
579
580
581 os << endl ;
582}
583
584
585////////////////////////////////////////////////////////////////////////////////
586/// Function called by RooMinimizer
587
588void RooFitResult::fillCorrMatrix(const std::vector<double>& globalCC, const TMatrixDSym& corrs, const TMatrixDSym& covs)
589{
590 // Sanity check
591 if (globalCC.empty() || corrs.GetNoElements() < 1 || covs.GetNoElements() < 1) {
592 coutI(Minimization) << "RooFitResult::fillCorrMatrix: number of floating parameters is zero, correlation matrix not filled" << endl ;
593 return ;
594 }
595
596 if (!_initPars) {
597 coutE(Minimization) << "RooFitResult::fillCorrMatrix: ERROR: list of initial parameters must be filled first" << endl ;
598 return ;
599 }
600
601 // Delete eventual prevous correlation data holders
602 if (_CM) delete _CM ;
603 if (_VM) delete _VM ;
604 if (_GC) delete _GC ;
605
606 // Build holding arrays for correlation coefficients
607 _CM = new TMatrixDSym(corrs) ;
608 _VM = new TMatrixDSym(covs) ;
609 _GC = new TVectorD(_CM->GetNcols()) ;
610 for(int i=0 ; i<_CM->GetNcols() ; i++) {
611 (*_GC)[i] = globalCC[i] ;
612 }
613 //fillLegacyCorrMatrix() ;
614}
615
616
617
618
619
620////////////////////////////////////////////////////////////////////////////////
621/// Sanity check
622
624{
625 if (!_CM) return ;
626
627 // Delete eventual prevous correlation data holders
628 if (_globalCorr) delete _globalCorr ;
630
631 // Build holding arrays for correlation coefficients
632 _globalCorr = new RooArgList("globalCorrelations") ;
633
634 for(RooAbsArg* arg : *_initPars) {
635 // Create global correlation value holder
636 TString gcName("GC[") ;
637 gcName.Append(arg->GetName()) ;
638 gcName.Append("]") ;
639 TString gcTitle(arg->GetTitle()) ;
640 gcTitle.Append(" Global Correlation") ;
641 _globalCorr->addOwned(*(new RooRealVar(gcName.Data(),gcTitle.Data(),0.))) ;
642
643 // Create array with correlation holders for this parameter
644 TString name("C[") ;
645 name.Append(arg->GetName()) ;
646 name.Append(",*]") ;
647 RooArgList* corrMatrixRow = new RooArgList(name.Data()) ;
648 _corrMatrix.Add(corrMatrixRow) ;
649 for(RooAbsArg* arg2 : *_initPars) {
650
651 TString cName("C[") ;
652 cName.Append(arg->GetName()) ;
653 cName.Append(",") ;
654 cName.Append(arg2->GetName()) ;
655 cName.Append("]") ;
656 TString cTitle("Correlation between ") ;
657 cTitle.Append(arg->GetName()) ;
658 cTitle.Append(" and ") ;
659 cTitle.Append(arg2->GetName()) ;
660 corrMatrixRow->addOwned(*(new RooRealVar(cName.Data(),cTitle.Data(),0.))) ;
661 }
662 }
663
664 for (unsigned int i = 0; i < (unsigned int)_CM->GetNcols() ; ++i) {
665
666 // Find the next global correlation slot to fill, skipping fixed parameters
667 auto& gcVal = static_cast<RooRealVar&>((*_globalCorr)[i]);
668 gcVal.setVal((*_GC)(i)) ; // WVE FIX THIS
669
670 // Fill a row of the correlation matrix
671 auto corrMatrixCol = static_cast<RooArgList const&>(*_corrMatrix.At(i));
672 for (unsigned int it = 0; it < (unsigned int)_CM->GetNcols() ; ++it) {
673 auto& cVal = static_cast<RooRealVar&>(corrMatrixCol[it]);
674 double value = (*_CM)(i,it) ;
675 cVal.setVal(value);
676 (*_CM)(i,it) = value;
677 }
678 }
679}
680
681
682
683
684
685////////////////////////////////////////////////////////////////////////////////
686/// Internal utility method to extract the correlation matrix and the
687/// global correlation coefficients from the MINUIT memory buffer and
688/// fill the internal arrays.
689
691{
692 // Sanity check
693 if (gMinuit->fNpar < 1) {
694 coutI(Minimization) << "RooFitResult::fillCorrMatrix: number of floating parameters is zero, correlation matrix not filled" << endl ;
695 return ;
696 }
697
698 if (!_initPars) {
699 coutE(Minimization) << "RooFitResult::fillCorrMatrix: ERROR: list of initial parameters must be filled first" << endl ;
700 return ;
701 }
702
703 // Delete eventual prevous correlation data holders
704 if (_CM) delete _CM ;
705 if (_VM) delete _VM ;
706 if (_GC) delete _GC ;
707
708 // Build holding arrays for correlation coefficients
709 _CM = new TMatrixDSym(_initPars->getSize()) ;
710 _VM = new TMatrixDSym(_initPars->getSize()) ;
711 _GC = new TVectorD(_initPars->getSize()) ;
712
713 // Extract correlation information for MINUIT (code taken from TMinuit::mnmatu() )
714
715 // WVE: This code directly manipulates minuit internal workspace,
716 // if TMinuit code changes this may need updating
717 Int_t ndex, i, j, m, n, it /* nparm,id,ix */ ;
718 Int_t ndi, ndj /*, iso, isw2, isw5*/;
719 for (i = 1; i <= gMinuit->fNpar; ++i) {
720 ndi = i*(i + 1) / 2;
721 for (j = 1; j <= gMinuit->fNpar; ++j) {
722 m = TMath::Max(i,j);
723 n = TMath::Min(i,j);
724 ndex = m*(m-1) / 2 + n;
725 ndj = j*(j + 1) / 2;
726 gMinuit->fMATUvline[j-1] = gMinuit->fVhmat[ndex-1] / TMath::Sqrt(TMath::Abs(gMinuit->fVhmat[ndi-1]*gMinuit->fVhmat[ndj-1]));
727 }
728
729 (*_GC)(i-1) = gMinuit->fGlobcc[i-1] ;
730
731 // Fill a row of the correlation matrix
732 for (it = 1; it <= gMinuit->fNpar ; ++it) {
733 (*_CM)(i-1,it-1) = gMinuit->fMATUvline[it-1] ;
734 }
735 }
736
737 for (int ii=0 ; ii<_finalPars->getSize() ; ii++) {
738 for (int jj=0 ; jj<_finalPars->getSize() ; jj++) {
739 (*_VM)(ii,jj) = (*_CM)(ii,jj) * ((RooRealVar*)_finalPars->at(ii))->getError() * ((RooRealVar*)_finalPars->at(jj))->getError() ;
740 }
741 }
742}
743
744////////////////////////////////////////////////////////////////////////////////
745
747{
748
749 // Delete eventual prevous correlation data holders
750 if (_CM)
751 delete _CM;
752 if (_VM)
753 delete _VM;
754 if (_GC)
755 delete _GC;
756
757 // Build holding arrays for correlation coefficients
760 _GC = new TVectorD(_initPars->getSize());
761
762 for (int ii = 0; ii < _finalPars->getSize(); ii++) {
763 (*_CM)(ii, ii) = 1;
764 (*_VM)(ii, ii) = ((RooRealVar *)_finalPars->at(ii))->getError() * ((RooRealVar *)_finalPars->at(ii))->getError();
765 (*_GC)(ii) = 0;
766 }
767}
768
769
770namespace {
771
772void isIdenticalErrMsg(std::string const& msgHead, const RooAbsReal* tv, const RooAbsReal* ov, bool verbose) {
773 if(!verbose) return;
774 std::cout << "RooFitResult::isIdentical: " << msgHead << " " << tv->GetName() << " differs in value:\t"
775 << tv->getVal() << " vs.\t" << ov->getVal()
776 << "\t(" << (tv->getVal()-ov->getVal())/ov->getVal() << ")" << std::endl;
777}
778
779void isErrorIdenticalErrMsg(std::string const& msgHead, const RooRealVar* tv, const RooRealVar* ov, bool verbose) {
780 if(!verbose) return;
781 std::cout << "RooFitResult::isIdentical: " << msgHead << " " << tv->GetName() << " differs in error:\t"
782 << tv->getError() << " vs.\t" << ov->getError()
783 << "\t(" << (tv->getError()-ov->getError())/ov->getError() << ")" << std::endl;
784}
785
786} // namespace
787
788
789////////////////////////////////////////////////////////////////////////////////
790/// Return true if this fit result is identical to other within tolerances, ignoring the correlation matrix.
791/// \param[in] other Fit result to test against.
792/// \param[in] tol **Relative** tolerance for parameters and NLL.
793/// \param[in] tolErr **Relative** tolerance for parameter errors.
794/// \param[in] verbose If this function will log to the standard output when comparisions fail.
795
796bool RooFitResult::isIdenticalNoCov(const RooFitResult& other, double tol, double tolErr, bool verbose) const
797{
798 bool ret = true;
799 auto deviation = [](const double left, const double right, double tolerance){
800 return right != 0. ? std::abs((left - right)/right) >= tolerance : std::abs(left) >= tolerance;
801 };
802
803 auto compare = [&](RooArgList const& pars, RooArgList const& otherpars, std::string const& prefix, bool isVerbose) {
804 bool out = true;
805
806 for (auto * tv : static_range_cast<const RooAbsReal*>(pars)) {
807 auto ov = static_cast<const RooAbsReal*>(otherpars.find(tv->GetName())) ;
808
809 // Check in the parameter is in the other fit result
810 if (!ov) {
811 if(verbose) cout << "RooFitResult::isIdentical: cannot find " << prefix << " " << tv->GetName() << " in reference" << endl ;
812 out = false;
813 }
814
815 // Compare parameter value
816 if (ov && deviation(tv->getVal(), ov->getVal(), tol)) {
817 isIdenticalErrMsg(prefix, tv, ov, isVerbose);
818 out = false;
819 }
820
821 // Compare parameter error if it's a RooRealVar
822 auto * rtv = dynamic_cast<RooRealVar const*>(tv);
823 auto * rov = dynamic_cast<RooRealVar const*>(ov);
824 if(rtv && rov) {
825 if (ov && deviation(rtv->getError(), rov->getError(), tolErr)) {
826 isErrorIdenticalErrMsg(prefix, rtv, rov, isVerbose);
827 out = false;
828 }
829 }
830 }
831
832 return out;
833 };
834
835 if (deviation(_minNLL, other._minNLL, tol)) {
836 if(verbose) std::cout << "RooFitResult::isIdentical: minimized value of -log(L) is different " << _minNLL << " vs. " << other._minNLL << std::endl;
837 ret = false;
838 }
839
840 ret &= compare(*_constPars, *other._constPars, "constant parameter", verbose);
841 ret &= compare(*_initPars, *other._initPars, "initial parameter", verbose);
842 ret &= compare(*_finalPars, *other._finalPars, "final parameter", verbose);
843
844 return ret;
845}
846
847
848////////////////////////////////////////////////////////////////////////////////
849/// Return true if this fit result is identical to other within tolerances.
850/// \param[in] other Fit result to test against.
851/// \param[in] tol **Relative** tolerance for parameters and NLL.
852/// \param[in] tolCorr **absolute** tolerance for correlation coefficients.
853/// \param[in] verbose If this function will log to the standard output when comparisions fail.
854///
855/// As the relative tolerance for the parameter errors, the default value of
856/// `1e-3` will be used.
857
858bool RooFitResult::isIdentical(const RooFitResult& other, double tol, double tolCorr, bool verbose) const
859{
860 bool ret = isIdenticalNoCov(other, tol, 1e-3 /* synced with default parameter*/, verbose);
861
862 auto deviationCorr = [tolCorr](const double left, const double right){
863 return std::abs(left - right) >= tolCorr;
864 };
865
866 // Only examine correlations for cases with >1 floating parameter
867 if (_finalPars->getSize()>1) {
868
870 other.fillLegacyCorrMatrix() ;
871
872 for (Int_t i=0 ; i<_globalCorr->getSize() ; i++) {
873 auto tv = static_cast<const RooAbsReal*>(_globalCorr->at(i));
874 auto ov = static_cast<const RooAbsReal*>(other._globalCorr->find(_globalCorr->at(i)->GetName())) ;
875 if (!ov) {
876 if(verbose) cout << "RooFitResult::isIdentical: cannot find global correlation coefficient " << tv->GetName() << " in reference" << endl ;
877 ret = false ;
878 }
879 if (ov && deviationCorr(tv->getVal(), ov->getVal())) {
880 isIdenticalErrMsg("global correlation coefficient", tv, ov, verbose);
881 ret = false ;
882 }
883 }
884
885 for (Int_t j=0 ; j<_corrMatrix.GetSize() ; j++) {
886 RooArgList* row = (RooArgList*) _corrMatrix.At(j) ;
887 RooArgList* orow = (RooArgList*) other._corrMatrix.At(j) ;
888 for (Int_t i=0 ; i<row->getSize() ; i++) {
889 auto tv = static_cast<const RooAbsReal*>(row->at(i));
890 auto ov = static_cast<const RooAbsReal*>(orow->find(tv->GetName())) ;
891 if (!ov) {
892 if(verbose) cout << "RooFitResult::isIdentical: cannot find correlation coefficient " << tv->GetName() << " in reference" << endl ;
893 ret = false ;
894 }
895 if (ov && deviationCorr(tv->getVal(), ov->getVal())) {
896 isIdenticalErrMsg("correlation coefficient", tv, ov, verbose);
897 ret = false ;
898 }
899 }
900 }
901 }
902
903 return ret ;
904}
905
906
907
908////////////////////////////////////////////////////////////////////////////////
909/// Import the results of the last fit performed by gMinuit, interpreting
910/// the fit parameters as the given varList of parameters.
911
913{
914 // Verify length of supplied varList
915 if (varList.getSize()>0 && varList.getSize()!=gMinuit->fNu) {
916 oocoutE(nullptr,InputArguments) << "RooFitResult::lastMinuitFit: ERROR: supplied variable list must be either empty " << endl
917 << " or match the number of variables of the last fit (" << gMinuit->fNu << ")" << endl ;
918 return nullptr;
919 }
920
921 // Verify that all members of varList are of type RooRealVar
922 for(RooAbsArg* arg : varList) {
923 if (!dynamic_cast<RooRealVar*>(arg)) {
924 oocoutE(nullptr,InputArguments) << "RooFitResult::lastMinuitFit: ERROR: variable '" << arg->GetName() << "' is not of type RooRealVar" << endl ;
925 return nullptr;
926 }
927 }
928
929 RooFitResult* r = new RooFitResult("lastMinuitFit","Last MINUIT fit") ;
930
931 // Extract names of fit parameters from MINUIT
932 // and construct corresponding RooRealVars
933 RooArgList constPars("constPars") ;
934 RooArgList floatPars("floatPars") ;
935
936 Int_t i ;
937 for (i = 1; i <= gMinuit->fNu; ++i) {
938 if (gMinuit->fNvarl[i-1] < 0) continue;
939 Int_t l = gMinuit->fNiofex[i-1];
940 TString varName(gMinuit->fCpnam[i-1]) ;
941 bool isConst(l==0) ;
942
943 double xlo = gMinuit->fAlim[i-1];
944 double xhi = gMinuit->fBlim[i-1];
945 double xerr = gMinuit->fWerr[l-1];
946 double xval = gMinuit->fU[i-1] ;
947
948 RooRealVar* var ;
949 if (varList.empty()) {
950
951 if ((xlo<xhi) && !isConst) {
952 var = new RooRealVar(varName,varName,xval,xlo,xhi) ;
953 } else {
954 var = new RooRealVar(varName,varName,xval) ;
955 }
956 var->setConstant(isConst) ;
957 } else {
958
959 var = (RooRealVar*) varList.at(i-1)->Clone() ;
960 var->setConstant(isConst) ;
961 var->setVal(xval) ;
962 if (xlo<xhi) {
963 var->setRange(xlo,xhi) ;
964 }
965 if (varName.CompareTo(var->GetName())) {
966 oocoutI(nullptr,Eval) << "RooFitResult::lastMinuitFit: fit parameter '" << varName
967 << "' stored in variable '" << var->GetName() << "'" << endl ;
968 }
969
970 }
971
972 if (isConst) {
973 constPars.addOwned(*var) ;
974 } else {
975 var->setError(xerr) ;
976 floatPars.addOwned(*var) ;
977 }
978 }
979
980 Int_t icode,npari,nparx ;
981 double fmin,edm,errdef ;
982 gMinuit->mnstat(fmin,edm,errdef,npari,nparx,icode) ;
983
984 r->setConstParList(constPars) ;
985 r->setInitParList(floatPars) ;
986 r->setFinalParList(floatPars) ;
987 r->setMinNLL(fmin) ;
988 r->setEDM(edm) ;
989 r->setCovQual(icode) ;
990 r->setStatus(gMinuit->fStatus) ;
991 r->fillCorrMatrix() ;
992
993 return r ;
994}
995
996
997
998////////////////////////////////////////////////////////////////////////////////
999/// Import the results of the last fit performed by gMinuit, interpreting
1000/// the fit parameters as the given varList of parameters.
1001
1003{
1004 // Verify that all members of varList are of type RooRealVar
1005 for(RooAbsArg * arg : paramList) {
1006 if (!dynamic_cast<RooRealVar *>(arg)) {
1007 oocoutE(nullptr, InputArguments) << "RooFitResult::lastMinuitFit: ERROR: variable '" << arg->GetName()
1008 << "' is not of type RooRealVar" << endl;
1009 return nullptr;
1010 }
1011 }
1012
1013 RooFitResult *r = new RooFitResult("lastMinuitFit", "Last MINUIT fit");
1014
1015 // Extract names of fit parameters from MINUIT
1016 // and construct corresponding RooRealVars
1017 RooArgList constPars("constPars");
1018 RooArgList floatPars("floatPars");
1019
1020 for(RooAbsArg* arg : paramList) {
1021 if (arg->isConstant()) {
1022 constPars.addClone(*arg);
1023 } else {
1024 floatPars.addClone(*arg);
1025 }
1026 }
1027
1028 r->setConstParList(constPars);
1029 r->setInitParList(floatPars);
1030 r->setFinalParList(floatPars);
1031 r->setMinNLL(0);
1032 r->setEDM(0);
1033 r->setCovQual(0);
1034 r->setStatus(0);
1035 r->fillPrefitCorrMatrix();
1036
1037 return r;
1038}
1039
1040////////////////////////////////////////////////////////////////////////////////
1041/// Store externally provided correlation matrix in this RooFitResult ;
1042
1044{
1045 // Delete any previous matrices
1046 if (_VM) {
1047 delete _VM ;
1048 }
1049 if (_CM) {
1050 delete _CM ;
1051 }
1052
1053 // Clone input covariance matrix ;
1054 _VM = (TMatrixDSym*) V.Clone() ;
1055
1056 // Now construct correlation matrix from it
1057 _CM = (TMatrixDSym*) _VM->Clone() ;
1058 for (Int_t i=0 ; i<_CM->GetNrows() ; i++) {
1059 for (Int_t j=0 ; j<_CM->GetNcols() ; j++) {
1060 if (i!=j) {
1061 (*_CM)(i,j) = (*_CM)(i,j) / sqrt((*_CM)(i,i)*(*_CM)(j,j)) ;
1062 }
1063 }
1064 }
1065 for (Int_t i=0 ; i<_CM->GetNrows() ; i++) {
1066 (*_CM)(i,i) = 1.0 ;
1067 }
1068
1069 _covQual = -1 ;
1070}
1071
1072
1073
1074////////////////////////////////////////////////////////////////////////////////
1075/// Return TH2D of correlation matrix
1076
1078{
1079 Int_t n = _CM->GetNcols() ;
1080
1081 TH2D* hh = new TH2D(name,name,n,0,n,n,0,n) ;
1082
1083 for (Int_t i = 0 ; i<n ; i++) {
1084 for (Int_t j = 0 ; j<n; j++) {
1085 hh->Fill(i+0.5,n-j-0.5,(*_CM)(i,j)) ;
1086 }
1087 hh->GetXaxis()->SetBinLabel(i+1,_finalPars->at(i)->GetName()) ;
1088 hh->GetYaxis()->SetBinLabel(n-i,_finalPars->at(i)->GetName()) ;
1089 }
1090 hh->SetMinimum(-1) ;
1091 hh->SetMaximum(+1) ;
1092
1093
1094 return hh ;
1095}
1096
1097
1098
1099
1100////////////////////////////////////////////////////////////////////////////////
1101/// Return covariance matrix
1102
1104{
1105 return *_VM ;
1106}
1107
1108
1109
1110
1111////////////////////////////////////////////////////////////////////////////////
1112/// Return a reduced covariance matrix (Note that Vred _is_ a simple sub-matrix of V,
1113/// row/columns are ordered to matched the convention given in input argument 'params'
1114
1116{
1117 const TMatrixDSym& V = covarianceMatrix() ;
1118
1119
1120 // Make sure that all given params were floating parameters in the represented fit
1121 RooArgList params2 ;
1122 for(RooAbsArg* arg : params) {
1123 if (_finalPars->find(arg->GetName())) {
1124 params2.add(*arg) ;
1125 } else {
1126 coutW(InputArguments) << "RooFitResult::reducedCovarianceMatrix(" << GetName() << ") WARNING input variable "
1127 << arg->GetName() << " was not a floating parameters in fit result and is ignored" << endl ;
1128 }
1129 }
1130
1131 // fix for bug ROOT-8044
1132 // use same order given bby vector params
1133 vector<int> indexMap(params2.getSize());
1134 for (int i=0 ; i<params2.getSize() ; i++) {
1135 indexMap[i] = _finalPars->index(params2[i].GetName());
1136 assert(indexMap[i] < V.GetNrows());
1137 }
1138
1139 TMatrixDSym Vred(indexMap.size());
1140 for (int i = 0; i < Vred.GetNrows(); ++i) {
1141 for (int j = 0; j < Vred.GetNcols(); ++j) {
1142 Vred(i,j) = V( indexMap[i], indexMap[j]);
1143 }
1144 }
1145 return Vred;
1146}
1147
1148
1149
1150////////////////////////////////////////////////////////////////////////////////
1151/// Return a reduced covariance matrix, which is calculated as
1152/// \f[
1153/// V_\mathrm{red} = \bar{V_{22}} = V_{11} - V_{12} \cdot V_{22}^{-1} \cdot V_{21},
1154/// \f]
1155/// where \f$ V_{11},V_{12},V_{21},V_{22} \f$ represent a block decomposition of the covariance matrix into observables that
1156/// are propagated (labeled by index '1') and that are not propagated (labeled by index '2'), and \f$ \bar{V_{22}} \f$
1157/// is the Shur complement of \f$ V_{22} \f$, calculated as shown above.
1158///
1159/// (Note that \f$ V_\mathrm{red} \f$ is *not* a simple sub-matrix of \f$ V \f$)
1160
1162{
1163 const TMatrixDSym& V = covarianceMatrix() ;
1164
1165 // Handle case where V==Vred here
1166 if (V.GetNcols()==params.getSize()) {
1167 return V ;
1168 }
1169
1170 double det = V.Determinant() ;
1171
1172 if (det<=0) {
1173 coutE(Eval) << "RooFitResult::conditionalCovarianceMatrix(" << GetName() << ") ERROR: covariance matrix is not positive definite (|V|="
1174 << det << ") cannot reduce it" << endl ;
1175 throw string("RooFitResult::conditionalCovarianceMatrix() ERROR, input covariance matrix is not positive definite") ;
1176 }
1177
1178 // Make sure that all given params were floating parameters in the represented fit
1179 RooArgList params2 ;
1180 for(RooAbsArg* arg : params) {
1181 if (_finalPars->find(arg->GetName())) {
1182 params2.add(*arg) ;
1183 } else {
1184 coutW(InputArguments) << "RooFitResult::conditionalCovarianceMatrix(" << GetName() << ") WARNING input variable "
1185 << arg->GetName() << " was not a floating parameters in fit result and is ignored" << endl ;
1186 }
1187 }
1188
1189 // Need to order params in vector in same order as in covariance matrix
1190 RooArgList params3 ;
1191 for(RooAbsArg* arg : *_finalPars) {
1192 if (params2.find(arg->GetName())) {
1193 params3.add(*arg) ;
1194 }
1195 }
1196
1197 // Find (subset) of parameters that are stored in the covariance matrix
1198 vector<int> map1, map2 ;
1199 for (int i=0 ; i<_finalPars->getSize() ; i++) {
1200 if (params3.find(_finalPars->at(i)->GetName())) {
1201 map1.push_back(i) ;
1202 } else {
1203 map2.push_back(i) ;
1204 }
1205 }
1206
1207 // Rearrange matrix in block form with 'params' first and 'others' last
1208 // (preserving relative order)
1209 TMatrixDSym S11, S22 ;
1210 TMatrixD S12, S21 ;
1211 RooMultiVarGaussian::blockDecompose(V,map1,map2,S11,S12,S21,S22) ;
1212
1213 // Constructed conditional matrix form -1
1214 // F(X1|X2) --> CovI --> S22bar = S11 - S12 S22 S21
1215
1216 // Do eigenvalue decomposition
1217 TMatrixD S22Inv(TMatrixD::kInverted,S22) ;
1218 TMatrixD S22bar = S11 - S12 * (S22Inv * S21) ;
1219
1220 // Convert explicitly to symmetric form
1221 TMatrixDSym Vred(S22bar.GetNcols()) ;
1222 for (int i=0 ; i<Vred.GetNcols() ; i++) {
1223 for (int j=i ; j<Vred.GetNcols() ; j++) {
1224 Vred(i,j) = (S22bar(i,j) + S22bar(j,i))/2 ;
1225 Vred(j,i) = Vred(i,j) ;
1226 }
1227 }
1228
1229 return Vred ;
1230}
1231
1232
1233
1234////////////////////////////////////////////////////////////////////////////////
1235/// Return correlation matrix ;
1236
1238{
1239 return *_CM ;
1240}
1241
1242
1243
1244////////////////////////////////////////////////////////////////////////////////
1245/// Return a p.d.f that represents the fit result as a multi-variate probability densisty
1246/// function on the floating fit parameters, including correlations
1247
1249{
1250 const TMatrixDSym& V = covarianceMatrix() ;
1251 double det = V.Determinant() ;
1252
1253 if (det<=0) {
1254 coutE(Eval) << "RooFitResult::createHessePdf(" << GetName() << ") ERROR: covariance matrix is not positive definite (|V|="
1255 << det << ") cannot construct p.d.f" << endl ;
1256 return 0 ;
1257 }
1258
1259 // Make sure that all given params were floating parameters in the represented fit
1260 RooArgList params2 ;
1261 for(RooAbsArg* arg : params) {
1262 if (_finalPars->find(arg->GetName())) {
1263 params2.add(*arg) ;
1264 } else {
1265 coutW(InputArguments) << "RooFitResult::createHessePdf(" << GetName() << ") WARNING input variable "
1266 << arg->GetName() << " was not a floating parameters in fit result and is ignored" << endl ;
1267 }
1268 }
1269
1270 // Need to order params in vector in same order as in covariance matrix
1271 RooArgList params3 ;
1272 for(RooAbsArg* arg : *_finalPars) {
1273 if (params2.find(arg->GetName())) {
1274 params3.add(*arg) ;
1275 }
1276 }
1277
1278
1279 // Handle special case of representing full covariance matrix here
1280 if (params3.getSize()==_finalPars->getSize()) {
1281
1282 RooArgList mu ;
1283 for (Int_t i=0 ; i<_finalPars->getSize() ; i++) {
1284 RooRealVar* parclone = (RooRealVar*) _finalPars->at(i)->Clone(Form("%s_centralvalue",_finalPars->at(i)->GetName())) ;
1285 parclone->setConstant(true) ;
1286 mu.add(*parclone) ;
1287 }
1288
1289 string name = Form("pdf_%s",GetName()) ;
1290 string title = Form("P.d.f of %s",GetTitle()) ;
1291
1292 // Create p.d.f.
1293 RooAbsPdf* mvg = new RooMultiVarGaussian(name.c_str(),title.c_str(),params3,mu,V) ;
1294 mvg->addOwnedComponents(mu) ;
1295 return mvg ;
1296 }
1297
1298 // -> ->
1299 // Handle case of conditional p.d.f. MVG(p1|p2) here
1300
1301 // Find (subset) of parameters that are stored in the covariance matrix
1302 vector<int> map1, map2 ;
1303 for (int i=0 ; i<_finalPars->getSize() ; i++) {
1304 if (params3.find(_finalPars->at(i)->GetName())) {
1305 map1.push_back(i) ;
1306 } else {
1307 map2.push_back(i) ;
1308 }
1309 }
1310
1311 // Rearrange matrix in block form with 'params' first and 'others' last
1312 // (preserving relative order)
1313 TMatrixDSym S11, S22 ;
1314 TMatrixD S12, S21 ;
1315 RooMultiVarGaussian::blockDecompose(V,map1,map2,S11,S12,S21,S22) ;
1316
1317 // Calculate offset vectors mu1 and mu2
1318 RooArgList mu1 ;
1319 for (UInt_t i=0 ; i<map1.size() ; i++) {
1320 RooRealVar* parclone = (RooRealVar*) _finalPars->at(map1[i])->Clone(Form("%s_centralvalue",_finalPars->at(map1[i])->GetName())) ;
1321 parclone->setConstant(true) ;
1322 mu1.add(*parclone) ;
1323 }
1324
1325 // Constructed conditional matrix form -1
1326 // F(X1|X2) --> CovI --> S22bar = S11 - S12 S22 S21
1327
1328 // Do eigenvalue decomposition
1329 TMatrixD S22Inv(TMatrixD::kInverted,S22) ;
1330 TMatrixD S22bar = S11 - S12 * (S22Inv * S21) ;
1331
1332 // Convert explicitly to symmetric form
1333 TMatrixDSym Vred(S22bar.GetNcols()) ;
1334 for (int i=0 ; i<Vred.GetNcols() ; i++) {
1335 for (int j=i ; j<Vred.GetNcols() ; j++) {
1336 Vred(i,j) = (S22bar(i,j) + S22bar(j,i))/2 ;
1337 Vred(j,i) = Vred(i,j) ;
1338 }
1339 }
1340 string name = Form("pdf_%s",GetName()) ;
1341 string title = Form("P.d.f of %s",GetTitle()) ;
1342
1343 // Create p.d.f.
1344 RooAbsPdf* ret = new RooMultiVarGaussian(name.c_str(),title.c_str(),params3,mu1,Vred) ;
1345 ret->addOwnedComponents(mu1) ;
1346 return ret ;
1347}
1348
1349
1350
1351////////////////////////////////////////////////////////////////////////////////
1352/// Change name of RooFitResult object
1353
1355{
1356 if (_dir) _dir->GetList()->Remove(this);
1358 if (_dir) _dir->GetList()->Add(this);
1359}
1360
1361
1362////////////////////////////////////////////////////////////////////////////////
1363/// Change name and title of RooFitResult object
1364
1365void RooFitResult::SetNameTitle(const char *name, const char* title)
1366{
1367 if (_dir) _dir->GetList()->Remove(this);
1368 TNamed::SetNameTitle(name,title) ;
1369 if (_dir) _dir->GetList()->Add(this);
1370}
1371
1372
1373////////////////////////////////////////////////////////////////////////////////
1374/// Print name of fit result
1375
1376void RooFitResult::printName(ostream& os) const
1377{
1378 os << GetName() ;
1379}
1380
1381
1382////////////////////////////////////////////////////////////////////////////////
1383/// Print title of fit result
1384
1385void RooFitResult::printTitle(ostream& os) const
1386{
1387 os << GetTitle() ;
1388}
1389
1390
1391////////////////////////////////////////////////////////////////////////////////
1392/// Print class name of fit result
1393
1394void RooFitResult::printClassName(ostream& os) const
1395{
1396 os << ClassName() ;
1397}
1398
1399
1400////////////////////////////////////////////////////////////////////////////////
1401/// Print arguments of fit result, i.e. the parameters of the fit
1402
1403void RooFitResult::printArgs(ostream& os) const
1404{
1405 os << "[constPars=" << *_constPars << ",floatPars=" << *_finalPars << "]" ;
1406}
1407
1408
1409
1410////////////////////////////////////////////////////////////////////////////////
1411/// Print the value of the fit result, i.e.g the status, minimized FCN, edm and covariance quality code
1412
1413void RooFitResult::printValue(ostream& os) const
1414{
1415 os << "(status=" << _status << ",FCNmin=" << _minNLL << ",EDM=" << _edm << ",covQual=" << _covQual << ")" ;
1416}
1417
1418
1419////////////////////////////////////////////////////////////////////////////////
1420/// Configure default contents to be printed
1421
1423{
1424 return kName|kClassName|kArgs|kValue ;
1425}
1426
1427
1428////////////////////////////////////////////////////////////////////////////////
1429/// Configure mapping of Print() arguments to RooPrintable print styles
1430
1432{
1433 if (!opt || strlen(opt)==0) {
1434 return kStandard ;
1435 }
1437}
1438
1439
1440////////////////////////////////////////////////////////////////////////////////
1441/// Stream an object of class RooFitResult.
1442
1444{
1445 if (R__b.IsReading()) {
1446 UInt_t R__s, R__c;
1447 Version_t R__v = R__b.ReadVersion(&R__s, &R__c);
1448 if (R__v>3) {
1449 R__b.ReadClassBuffer(RooFitResult::Class(),this,R__v,R__s,R__c);
1452 } else {
1453 // backward compatibitily streaming
1454 TNamed::Streamer(R__b);
1457 R__b >> _status;
1458 R__b >> _covQual;
1459 R__b >> _numBadNLL;
1460 R__b >> _minNLL;
1461 R__b >> _edm;
1462 R__b >> _constPars;
1463 R__b >> _initPars;
1464 R__b >> _finalPars;
1465 R__b >> _globalCorr;
1466 _corrMatrix.Streamer(R__b);
1467 R__b.CheckByteCount(R__s, R__c, RooFitResult::IsA());
1468
1469 // Now fill new-style covariance and correlation matrix information
1470 // from legacy form
1471 _CM = new TMatrixDSym(_finalPars->getSize()) ;
1472 _VM = new TMatrixDSym(_CM->GetNcols()) ;
1473 _GC = new TVectorD(_CM->GetNcols()) ;
1474
1475 for (unsigned int i = 0; i < (unsigned int)_CM->GetNcols() ; ++i) {
1476
1477 // Find the next global correlation slot to fill, skipping fixed parameters
1478 auto& gcVal = static_cast<RooRealVar&>((*_globalCorr)[i]);
1479 (*_GC)(i) = gcVal.getVal() ;
1480
1481 // Fill a row of the correlation matrix
1482 auto corrMatrixCol = static_cast<RooArgList const&>(*_corrMatrix.At(i));
1483 for (unsigned int it = 0; it < (unsigned int)_CM->GetNcols() ; ++it) {
1484 auto& cVal = static_cast<RooRealVar&>(corrMatrixCol[it]);
1485 double value = cVal.getVal() ;
1486 (*_CM)(it,i) = value ;
1487 (*_CM)(i,it) = value;
1488 (*_VM)(it,i) = value*((RooRealVar*)_finalPars->at(i))->getError()*((RooRealVar*)_finalPars->at(it))->getError() ;
1489 (*_VM)(i,it) = (*_VM)(it,i) ;
1490 }
1491 }
1492 }
1493
1494 } else {
1496 }
1497}
1498
#define s1(x)
Definition: RSha256.hxx:91
#define e(i)
Definition: RSha256.hxx:103
#define coutI(a)
Definition: RooMsgService.h:34
#define coutW(a)
Definition: RooMsgService.h:36
#define oocoutE(o, a)
Definition: RooMsgService.h:52
#define oocoutI(o, a)
Definition: RooMsgService.h:49
#define coutE(a)
Definition: RooMsgService.h:37
short Version_t
Definition: RtypesCore.h:65
const char Option_t
Definition: RtypesCore.h:66
#define ClassImp(name)
Definition: Rtypes.h:375
@ kRed
Definition: Rtypes.h:66
@ kBlack
Definition: Rtypes.h:65
@ kBlue
Definition: Rtypes.h:66
@ kDashed
Definition: TAttLine.h:48
static void indent(ostringstream &buf, int indent_level)
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 r
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 value
Option_t Option_t TPoint TPoint const char x2
Option_t Option_t TPoint TPoint const char x1
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 g
char name[80]
Definition: TGX11.cxx:110
TMatrixTSym< Double_t > TMatrixDSym
TMatrixT< Float_t > TMatrix
Definition: TMatrix.h:24
R__EXTERN TMinuit * gMinuit
Definition: TMinuit.h:271
char * Form(const char *fmt,...)
Formats a string in a circular formatting buffer.
Definition: TString.cxx:2468
TVectorT< Double_t > TVectorD
Definition: TVectorDfwd.h:23
RooAbsArg is the common abstract base class for objects that represent a value and a "shape" in RooFi...
Definition: RooAbsArg.h:72
static void ioStreamerPass2Finalize()
Method called by workspace container to finalize schema evolution issues that cannot be handled in a ...
Definition: RooAbsArg.cxx:2419
bool addOwnedComponents(const RooAbsCollection &comps)
Take ownership of the contents of 'comps'.
Definition: RooAbsArg.cxx:2211
TObject * Clone(const char *newname=nullptr) const override
Make a clone of an object using the Streamer facility.
Definition: RooAbsArg.h:84
bool empty() const
RooAbsCollection * snapshot(bool deepCopy=true) const
Take a snap shot of current collection contents.
Int_t getSize() const
Return the number of elements in the collection.
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...
virtual bool addOwned(RooAbsArg &var, bool silent=false)
Add an argument and transfer the ownership to the collection.
virtual RooAbsArg * addClone(const RooAbsArg &var, bool silent=false)
Add a clone of the specified argument to list.
RooAbsArg * find(const char *name) const
Find object with given name in list.
void setConstant(bool value=true)
RooAbsReal is the common abstract base class for objects that represent a real value and implements f...
Definition: RooAbsReal.h:60
double getVal(const RooArgSet *normalisationSet=nullptr) const
Evaluate object.
Definition: RooAbsReal.h:104
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:56
RooDirItem is a utility base class for RooFit objects that are to be attached to ROOT directories.
Definition: RooDirItem.h:22
virtual void Streamer(TBuffer &)
void removeFromDir(TObject *obj)
Remove object from directory it was added to.
Definition: RooDirItem.cxx:39
TDirectory * _dir
! Associated directory
Definition: RooDirItem.h:33
void appendToDir(TObject *obj, bool forceMemoryResident=false)
Append object to directory.
Definition: RooDirItem.cxx:51
A RooEllipse is a two-dimensional ellipse that can be used to represent an error contour.
Definition: RooEllipse.h:22
RooFitResult is a container class to hold the input and output of a PDF fit to a dataset.
Definition: RooFitResult.h:40
TMatrixDSym conditionalCovarianceMatrix(const RooArgList &params) const
Return a reduced covariance matrix, which is calculated as.
void fillCorrMatrix()
Internal utility method to extract the correlation matrix and the global correlation coefficients fro...
double correlation(const RooAbsArg &par1, const RooAbsArg &par2) const
Return correlation between par1 and par2.
Definition: RooFitResult.h:118
~RooFitResult() override
Destructor.
TList _corrMatrix
! Correlation matrix (list of RooArgLists)
Definition: RooFitResult.h:192
std::vector< std::pair< std::string, int > > _statusHistory
History of status codes.
Definition: RooFitResult.h:201
const TMatrixDSym & covarianceMatrix() const
Return covariance matrix.
TMatrixDSym * _CM
Correlation matrix.
Definition: RooFitResult.h:197
bool isIdentical(const RooFitResult &other, double tol=1e-6, double tolCorr=1e-4, bool verbose=true) const
Return true if this fit result is identical to other within tolerances.
void setConstParList(const RooArgList &list)
Fill the list of constant parameters.
Int_t statusCodeHistory(UInt_t icycle) const
double _minNLL
NLL at minimum.
Definition: RooFitResult.h:185
Int_t _numBadNLL
Number calls with bad (zero,negative) likelihood.
Definition: RooFitResult.h:184
TMatrixDSym * _VM
Covariance matrix.
Definition: RooFitResult.h:198
void SetNameTitle(const char *name, const char *title) override
Change name and title of RooFitResult object.
Int_t defaultPrintContents(Option_t *opt) const override
Configure default contents to be printed.
void printClassName(std::ostream &os) const override
Print class name of fit result.
RooArgList * _initPars
List of floating parameters with initial values.
Definition: RooFitResult.h:188
TClass * IsA() const override
Definition: RooFitResult.h:203
RooFitResult(const char *name=nullptr, const char *title=nullptr)
Constructor with name and title.
TMatrixDSym reducedCovarianceMatrix(const RooArgList &params) const
Return a reduced covariance matrix (Note that Vred is a simple sub-matrix of V, row/columns are order...
Int_t _covQual
MINUIT quality code of covariance matrix.
Definition: RooFitResult.h:183
void fillPrefitCorrMatrix()
RooArgList * _constPars
List of constant parameters.
Definition: RooFitResult.h:187
RooArgList * _globalCorr
! List of global correlation coefficients
Definition: RooFitResult.h:191
const RooArgList & randomizePars() const
Generate random perturbations of the final parameters using the covariance matrix.
const RooArgList * globalCorr()
Return the list of all global correlations.
static RooFitResult * prefitResult(const RooArgList &paramList)
Import the results of the last fit performed by gMinuit, interpreting the fit parameters as the given...
void setCovarianceMatrix(TMatrixDSym &V)
Store externally provided correlation matrix in this RooFitResult ;.
RooArgList * _finalPars
List of floating parameters with final values.
Definition: RooFitResult.h:189
double edm() const
Return estimated distance to minimum.
Definition: RooFitResult.h:95
const RooArgList & constPars() const
Return list of constant parameters.
Definition: RooFitResult.h:103
Int_t _status
MINUIT status code.
Definition: RooFitResult.h:182
bool isIdenticalNoCov(const RooFitResult &other, double tol=1e-6, double tolErr=1e-3, bool verbose=true) const
Return true if this fit result is identical to other within tolerances, ignoring the correlation matr...
void SetName(const char *name) override
Change name of RooFitResult object.
const char * statusLabelHistory(UInt_t icycle) const
static TClass * Class()
void printValue(std::ostream &os) const override
Print the value of the fit result, i.e.g the status, minimized FCN, edm and covariance quality code.
void printTitle(std::ostream &os) const override
Print title of fit result.
TH2 * correlationHist(const char *name="correlation_matrix") const
Return TH2D of correlation matrix.
void fillLegacyCorrMatrix() const
Sanity check.
void setInitParList(const RooArgList &list)
Fill the list of initial values of the floating parameters.
RooPlot * plotOn(RooPlot *frame, const RooAbsArg &par1, const RooAbsArg &par2, const char *options="ME") const
Add objects to a 2D plot.
Definition: RooFitResult.h:144
void Streamer(TBuffer &) override
Stream an object of class RooFitResult.
StyleOption defaultPrintStyle(Option_t *opt) const override
Configure mapping of Print() arguments to RooPrintable print styles.
double covariance(Int_t row, Int_t col) const
Return the covariance matrix element addressed with numeric indices.
TMatrixF * _Lt
! triangular matrix used for generate random perturbations
Definition: RooFitResult.h:195
const RooArgList & floatParsFinal() const
Return list of floating parameters after fit.
Definition: RooFitResult.h:111
RooArgList * _randomPars
! List of floating parameters with most recent random perturbation applied
Definition: RooFitResult.h:194
static RooFitResult * lastMinuitFit(const RooArgList &varList=RooArgList())
Import the results of the last fit performed by gMinuit, interpreting the fit parameters as the given...
void printMultiline(std::ostream &os, Int_t contents, bool verbose=false, TString indent="") const override
Print fit result to stream 'os'.
TVectorD * _GC
Global correlation coefficients.
Definition: RooFitResult.h:199
RooAbsPdf * createHessePdf(const RooArgSet &params) const
Return a p.d.f that represents the fit result as a multi-variate probability densisty function on the...
double _edm
Estimated distance to minimum.
Definition: RooFitResult.h:186
void setFinalParList(const RooArgList &list)
Fill the list of final values of the floating parameters.
void printArgs(std::ostream &os) const override
Print arguments of fit result, i.e. the parameters of the fit.
void printName(std::ostream &os) const override
Print name of fit result.
const TMatrixDSym & correlationMatrix() const
Return correlation matrix ;.
Multivariate Gaussian p.d.f.
static void blockDecompose(const TMatrixD &input, const std::vector< int > &map1, const std::vector< int > &map2, TMatrixDSym &S11, TMatrixD &S12, TMatrixD &S21, TMatrixDSym &S22)
Block decomposition of covI according to given maps of observables.
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:383
void addPlotable(RooPlotable *plotable, Option_t *drawOptions="", bool invisible=false, bool refreshNorm=false)
Add the specified plotable object to our plot.
Definition: RooPlot.cxx:531
RooPlotable is a 'mix-in' base class that define the standard RooFit plotting and printing methods.
Definition: RooPrintable.h:25
virtual StyleOption defaultPrintStyle(Option_t *opt) const
virtual void Streamer(TBuffer &)
static double gaussian(TRandom *generator=randomGenerator())
Return a Gaussian random variable with mean 0 and variance 1.
Definition: RooRandom.cxx:108
RooRealVar represents a variable that can be changed from the outside.
Definition: RooRealVar.h:40
void setVal(double value) override
Set value of variable to 'value'.
Definition: RooRealVar.cxx:254
void setError(double value)
Definition: RooRealVar.h:64
double getError() const
Definition: RooRealVar.h:62
void setRange(const char *name, double min, double max)
Set a fit or plotting range.
Definition: RooRealVar.cxx:525
virtual void SetLineStyle(Style_t lstyle)
Set the line style.
Definition: TAttLine.h:42
virtual void SetLineWidth(Width_t lwidth)
Set the line width.
Definition: TAttLine.h:43
virtual void SetLineColor(Color_t lcolor)
Set the line color.
Definition: TAttLine.h:40
virtual void SetMarkerColor(Color_t mcolor=1)
Set the marker color.
Definition: TAttMarker.h:38
virtual void SetBinLabel(Int_t bin, const char *label)
Set label for bin.
Definition: TAxis.cxx:851
Create a Box.
Definition: TBox.h:22
Buffer base class used for serializing objects.
Definition: TBuffer.h:43
virtual Version_t ReadVersion(UInt_t *start=nullptr, UInt_t *bcnt=nullptr, const TClass *cl=nullptr)=0
virtual Int_t CheckByteCount(UInt_t startpos, UInt_t bcnt, const TClass *clss)=0
virtual Int_t ReadClassBuffer(const TClass *cl, void *pointer, const TClass *onfile_class=nullptr)=0
Bool_t IsReading() const
Definition: TBuffer.h:86
virtual Int_t WriteClassBuffer(const TClass *cl, void *pointer)=0
virtual void RemoveAll(TCollection *col)
Remove all objects in collection col from this collection.
virtual void SetOwner(Bool_t enable=kTRUE)
Set whether this collection is the owner (enable==true) of its content.
virtual Int_t GetSize() const
Return the capacity of the collection, i.e.
Definition: TCollection.h:184
virtual TList * GetList() const
Definition: TDirectory.h:222
The axis painter class.
Definition: TGaxis.h:23
TAxis * GetXaxis()
Definition: TH1.h:322
virtual void SetMaximum(Double_t maximum=-1111)
Definition: TH1.h:400
TAxis * GetYaxis()
Definition: TH1.h:323
virtual void SetMinimum(Double_t minimum=-1111)
Definition: TH1.h:401
2-D histogram with a double per channel (see TH1 documentation)}
Definition: TH2.h:300
Service class for 2-D histogram classes.
Definition: TH2.h:30
Int_t Fill(Double_t) override
Invalid Fill method.
Definition: TH2.cxx:347
Use the TLine constructor to create a simple line.
Definition: TLine.h:22
void Streamer(TBuffer &) override
Stream all objects in the collection to or from the I/O buffer.
Definition: TList.cxx:1191
void Add(TObject *obj) override
Definition: TList.h:81
TObject * Remove(TObject *obj) override
Remove object from the list.
Definition: TList.cxx:822
void Delete(Option_t *option="") override
Remove all objects from the list AND delete all heap based objects.
Definition: TList.cxx:470
TObject * At(Int_t idx) const override
Returns the object at position idx. Returns 0 if idx is out of range.
Definition: TList.cxx:357
Manages Markers.
Definition: TMarker.h:22
Int_t GetNrows() const
Definition: TMatrixTBase.h:123
Int_t GetNoElements() const
Definition: TMatrixTBase.h:127
Int_t GetNcols() const
Definition: TMatrixTBase.h:126
Double_t Determinant() const override
TMatrixT.
Definition: TMatrixT.h:39
@ kTransposed
Definition: TMatrixT.h:58
Double_t * fU
Definition: TMinuit.h:68
Int_t fNu
Definition: TMinuit.h:130
Double_t * fGlobcc
Definition: TMinuit.h:74
TString * fCpnam
Character to be plotted at the X,Y contour positions.
Definition: TMinuit.h:165
Int_t * fNvarl
Definition: TMinuit.h:126
Double_t * fMATUvline
Definition: TMinuit.h:107
Double_t * fBlim
Definition: TMinuit.h:70
Double_t * fWerr
Definition: TMinuit.h:73
Int_t fStatus
Definition: TMinuit.h:154
Int_t fNpar
Definition: TMinuit.h:41
Double_t * fAlim
Definition: TMinuit.h:69
Int_t * fNiofex
Definition: TMinuit.h:127
Double_t * fVhmat
Definition: TMinuit.h:89
virtual void mnstat(Double_t &fmin, Double_t &fedm, Double_t &errdef, Int_t &npari, Int_t &nparx, Int_t &istat)
Returns concerning the current status of the minimization.
Definition: TMinuit.cxx:7637
The TNamed class is the base class for all named ROOT classes.
Definition: TNamed.h:29
const char * GetName() const override
Returns name of object.
Definition: TNamed.h:47
void Streamer(TBuffer &) override
Stream an object of class TObject.
const char * GetTitle() const override
Returns title of object.
Definition: TNamed.h:48
virtual void SetName(const char *name)
Set the name of the TNamed.
Definition: TNamed.cxx:140
virtual void SetNameTitle(const char *name, const char *title)
Set all the TNamed parameters (name and title).
Definition: TNamed.cxx:154
virtual TObject * Clone(const char *newname="") const
Make a clone of an object using the Streamer facility.
Definition: TObject.cxx:223
virtual const char * ClassName() const
Returns name of class to which the object belongs.
Definition: TObject.cxx:207
Basic string class.
Definition: TString.h:136
int CompareTo(const char *cs, ECaseCompare cmp=kExact) const
Compare a string to char *cs2.
Definition: TString.cxx:451
const char * Data() const
Definition: TString.h:369
void ToUpper()
Change string to upper case.
Definition: TString.cxx:1184
TString & Append(const char *cs)
Definition: TString.h:565
Bool_t Contains(const char *pat, ECaseCompare cmp=kExact) const
Definition: TString.h:625
TVectorT.
Definition: TVectorT.h:27
TLine * line
void box(Int_t pat, Double_t x1, Double_t y1, Double_t x2, Double_t y2)
Definition: fillpatterns.C:1
RVec< PromoteType< T > > abs(const RVec< T > &v)
Definition: RVec.hxx:1778
const Int_t n
Definition: legend1.C:16
for(Int_t i=0;i< n;i++)
Definition: legend1.C:18
VecExpr< UnaryOp< Sqrt< T >, VecExpr< A, T, D >, T >, T, D > sqrt(const VecExpr< A, T, D > &rhs)
@ Minimization
Definition: RooGlobalFunc.h:62
@ InputArguments
Definition: RooGlobalFunc.h:63
static constexpr double L
Short_t Max(Short_t a, Short_t b)
Returns the largest of a and b.
Definition: TMathBase.h:250
Double_t Sqrt(Double_t x)
Returns the square root of x.
Definition: TMath.h:660
Short_t Min(Short_t a, Short_t b)
Returns the smallest of a and b.
Definition: TMathBase.h:198
Short_t Abs(Short_t d)
Returns the absolute value of parameter Short_t d.
Definition: TMathBase.h:123
TMarker m
Definition: textangle.C:8
TLine l
Definition: textangle.C:4