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