Logo ROOT  
Reference Guide
TFitResult.cxx
Go to the documentation of this file.
1 // @(#)root/hist:$Id$
2 // Author: David Gonzalez Maline 12/11/09
3 
4 /*************************************************************************
5  * Copyright (C) 1995-2000, Rene Brun and Fons Rademakers. *
6  * All rights reserved. *
7  * *
8  * For the licensing terms see $ROOTSYS/LICENSE. *
9  * For the list of contributors see $ROOTSYS/README/CREDITS. *
10  *************************************************************************/
11 
12 #include "TFitResult.h"
13 #include "Math/WrappedMultiTF1.h"
14 #include "TGraph.h"
15 
16 
17 #include <iostream>
18 
19 /** \class TFitResult
20  \ingroup Hist
21 Extends the ROOT::Fit::Result class with a TNamed inheritance
22 providing easy possibility for I/O
23 */
24 
26 
27 ////////////////////////////////////////////////////////////////////////////////
28 /// Constructor from a ROOT::Fit::FitResult
29 /// copy the contained TF1 pointer function if it is
30 
32  TNamed("TFitResult","TFitResult"),
33  ROOT::Fit::FitResult(f)
34 {
35  ROOT::Math::WrappedMultiTF1 * wfunc = dynamic_cast<ROOT::Math::WrappedMultiTF1 *>(ModelFunction().get() );
36  if (wfunc) wfunc->SetAndCopyFunction();
37 }
38 
39 
40 ////////////////////////////////////////////////////////////////////////////////
41 /// Print result of the fit, by default chi2, parameter values and errors.
42 /// if option "V" is given print also error matrix and correlation
43 
44 void TFitResult::Print(Option_t *option) const
45 {
46  TString opt(option);
47  opt.ToUpper();
48  bool doCovMat = opt.Contains("V");
49  ROOT::Fit::FitResult::Print( std::cout, doCovMat);
50 }
51 
52 ////////////////////////////////////////////////////////////////////////////////
53 /// Return the covariance matrix from fit
54 ///
55 /// The matrix is a symmetric matrix with a size N equal to
56 /// the total number of parameters considered in the fit including the fixed ones
57 /// The matrix row and columns corresponding to the fixed parameters will contain only zero's
58 
60 {
61  if (CovMatrixStatus() == 0) {
62  Warning("GetCovarianceMatrix","covariance matrix is not available");
63  return TMatrixDSym();
64  }
65  TMatrixDSym mat(NPar());
66  ROOT::Fit::FitResult::GetCovarianceMatrix<TMatrixDSym>(mat);
67  return mat;
68 }
69 
70 ////////////////////////////////////////////////////////////////////////////////
71 /// Return the correlation matrix from fit.
72 ///
73 /// The matrix is a symmetric matrix with a size N equal to
74 /// the total number of parameters considered in the fit including the fixed ones
75 /// The matrix row and columns corresponding to the fixed parameters will contain only zero's
76 
78 {
79  if (CovMatrixStatus() == 0) {
80  Warning("GetCorrelationMatrix","correlation matrix is not available");
81  return TMatrixDSym();
82  }
83  TMatrixDSym mat(NPar());
84  ROOT::Fit::FitResult::GetCorrelationMatrix<TMatrixDSym>(mat);
85  return mat;
86 }
87 
88 ////////////////////////////////////////////////////////////////////////////////
89 /// Scan parameter ipar between value of xmin and xmax
90 /// A graph must be given which will be on return filled with the scan resul
91 /// If the graph size is zero, a default size n = 40 will be used
92 
93 bool TFitResult::Scan(unsigned int ipar, TGraph *gr, double xmin, double xmax)
94 {
95  if (!gr)
96  return false;
97 
98  unsigned int npoints = gr->GetN();
99  if (npoints == 0) {
100  npoints = 40;
101  gr->Set(npoints);
102  }
103  bool ret = ROOT::Fit::FitResult::Scan(ipar, npoints, gr->GetX(), gr->GetY(), xmin, xmax);
104  if ((int)npoints < gr->GetN())
105  gr->Set(npoints);
106  return ret;
107 }
108 
109 ////////////////////////////////////////////////////////////////////////////////
110 /// Create a 2D contour around the minimum for the parameter ipar and jpar
111 /// if a minimum does not exist or is invalid it will return false
112 /// on exit a TGraph is filled with the contour points
113 /// the number of contour points is determined by the size of the TGraph.
114 /// if the size is zero a default number of points = 20 is used
115 /// pass optionally the confidence level, default is 0.683
116 /// it is assumed that ErrorDef() defines the right error definition
117 /// (i.e 1 sigma error for one parameter). If not the confidence level are scaled to new level
118 
119 bool TFitResult::Contour(unsigned int ipar, unsigned int jpar, TGraph *gr, double confLevel)
120 {
121  if (!gr)
122  return false;
123 
124  unsigned int npoints = gr->GetN();
125  if (npoints == 0) {
126  npoints = 40;
127  gr->Set(npoints);
128  }
129  bool ret = ROOT::Fit::FitResult::Contour(ipar, jpar, npoints, gr->GetX(), gr->GetY(), confLevel);
130  if ((int)npoints < gr->GetN())
131  gr->Set(npoints);
132 
133  return ret;
134 }
135 
136 ////////////////////////////////////////////////////////////////////////////////
137 /// Print the TFitResult.
138 
139 std::string cling::printValue(const TFitResult* val)
140 {
141  std::stringstream outs;
142  val->ROOT::Fit::FitResult::Print(outs, false /*doCovMat*/);
143  return outs.str();
144 }
ROOT::Math::WrappedMultiTF1Templ
Class to Wrap a ROOT Function class (like TF1) in a IParamMultiFunction interface of multi-dimensions...
Definition: WrappedMultiTF1.h:57
ROOT::Fit::FitResult::Print
void Print(std::ostream &os, bool covmat=false) const
print the result and optionaly covariance matrix and correlations
Definition: FitResult.cxx:438
TFitResult::GetCovarianceMatrix
TMatrixDSym GetCovarianceMatrix() const
Return the covariance matrix from fit.
Definition: TFitResult.cxx:59
HFit::Fit
TFitResultPtr Fit(FitObject *h1, TF1 *f1, Foption_t &option, const ROOT::Math::MinimizerOptions &moption, const char *goption, ROOT::Fit::DataRange &range)
Definition: HFitImpl.cxx:133
TFitResult.h
f
#define f(i)
Definition: RSha256.hxx:122
ROOT::Fit::FitResult::Contour
bool Contour(unsigned int ipar, unsigned int jpar, unsigned int &npoints, double *pntsx, double *pntsy, double confLevel=0.683)
create contour of two parameters around the minimum pass as option confidence level: default is a val...
Definition: FitResult.cxx:709
ClassImp
#define ClassImp(name)
Definition: Rtypes.h:364
TGraph.h
xmax
float xmax
Definition: THbookFile.cxx:95
TFitResult
Definition: TFitResult.h:34
TString::Contains
Bool_t Contains(const char *pat, ECaseCompare cmp=kExact) const
Definition: TString.h:624
TMatrixTSym
Definition: TMatrixDSymfwd.h:22
ROOT::Fit::FitResult::ModelFunction
std::shared_ptr< IModelFunction > ModelFunction()
Return pointer non const pointer to model (fit) function with fitted parameter values.
Definition: FitResult.h:343
TString
Definition: TString.h:136
TFitResult::GetCorrelationMatrix
TMatrixDSym GetCorrelationMatrix() const
Return the correlation matrix from fit.
Definition: TFitResult.cxx:77
ROOT::Fit::FitResult::NPar
unsigned int NPar() const
total number of parameters (abbreviation)
Definition: FitResult.h:136
TString::ToUpper
void ToUpper()
Change string to upper case.
Definition: TString.cxx:1138
ROOT::Fit::FitResult::CovMatrixStatus
int CovMatrixStatus() const
covariance matrix status code using Minuit convention : =0 not calculated, =1 approximated,...
Definition: FitResult.h:147
ROOT::Math::WrappedMultiTF1Templ::SetAndCopyFunction
void SetAndCopyFunction(const TF1 *f=0)
method to set a new function pointer and copy it inside.
Definition: WrappedMultiTF1.h:340
TGraph::GetX
Double_t * GetX() const
Definition: TGraph.h:130
Option_t
const typedef char Option_t
Definition: RtypesCore.h:66
ROOT::Fit::FitResult::Scan
bool Scan(unsigned int ipar, unsigned int &npoints, double *pntsx, double *pntsy, double xmin=0, double xmax=0)
scan likelihood value of parameter and fill the given graph.
Definition: FitResult.cxx:688
xmin
float xmin
Definition: THbookFile.cxx:95
gr
TGraphErrors * gr
Definition: legend1.C:25
TNamed
Definition: TNamed.h:29
ROOT::Fit::FitResult
class containg the result of the fit and all the related information (fitted parameter values,...
Definition: FitResult.h:52
TFitResult::Contour
bool Contour(unsigned int ipar, unsigned int jpar, TGraph *gr, double confLevel=0.683)
Create a 2D contour around the minimum for the parameter ipar and jpar if a minimum does not exist or...
Definition: TFitResult.cxx:119
TGraph::GetY
Double_t * GetY() const
Definition: TGraph.h:131
TFitResult::TFitResult
TFitResult(int status=0)
Definition: TFitResult.h:39
TObject::Warning
virtual void Warning(const char *method, const char *msgfmt,...) const
Issue warning message.
Definition: TObject.cxx:876
TFitResult::Print
virtual void Print(Option_t *option="") const
Print result of the fit, by default chi2, parameter values and errors.
Definition: TFitResult.cxx:44
TGraph::Set
virtual void Set(Int_t n)
Set number of points in the graph Existing coordinates are preserved New coordinates above fNpoints a...
Definition: TGraph.cxx:2204
TGraph
Definition: TGraph.h:41
WrappedMultiTF1.h
TGraph::GetN
Int_t GetN() const
Definition: TGraph.h:123
TFitResult::Scan
bool Scan(unsigned int ipar, TGraph *gr, double xmin=0, double xmax=0)
Scan parameter ipar between value of xmin and xmax A graph must be given which will be on return fill...
Definition: TFitResult.cxx:93
TMatrixDSym
TMatrixTSym< Double_t > TMatrixDSym
Definition: TMatrixDSymfwd.h:22
ROOT
VSD Structures.
Definition: StringConv.hxx:21