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