Logo ROOT  
Reference Guide
HFitImpl.cxx
Go to the documentation of this file.
1 // new HFit function
2 //______________________________________________________________________________
3 
4 
5 #include "TH1.h"
6 #include "TF1.h"
7 #include "TF2.h"
8 #include "TF3.h"
9 #include "TError.h"
10 #include "TGraph.h"
11 #include "TMultiGraph.h"
12 #include "TGraph2D.h"
13 #include "THnBase.h"
14 
15 #include "Fit/Fitter.h"
16 #include "Fit/FitConfig.h"
17 #include "Fit/BinData.h"
18 #include "Fit/UnBinData.h"
19 #include "Fit/Chi2FCN.h"
21 #include "HFitInterface.h"
22 #include "Math/MinimizerOptions.h"
23 #include "Math/Minimizer.h"
24 
25 #include "Math/WrappedTF1.h"
26 #include "Math/WrappedMultiTF1.h"
27 
28 #include "TList.h"
29 #include "TMath.h"
30 #include "TROOT.h"
31 
32 #include "TVirtualPad.h" // for gPad
33 
34 #include "TBackCompFitter.h"
35 #include "TFitResultPtr.h"
36 #include "TFitResult.h"
37 
38 #include <stdlib.h>
39 #include <cmath>
40 #include <memory>
41 #include <limits>
42 
43 //#define DEBUG
44 
45 // utility functions used in TH1::Fit
46 
47 namespace HFit {
48 
49 
50  int GetDimension(const TH1 * h1) { return h1->GetDimension(); }
51  int GetDimension(const TGraph * ) { return 1; }
52  int GetDimension(const TMultiGraph * ) { return 1; }
53  int GetDimension(const TGraph2D * ) { return 2; }
54  int GetDimension(const THnBase * s1) { return s1->GetNdimensions(); }
55 
56  int CheckFitFunction(const TF1 * f1, int hdim);
57 
58 
59  void GetFunctionRange(const TF1 & f1, ROOT::Fit::DataRange & range);
60 
61  void FitOptionsMake(const char *option, Foption_t &fitOption);
62 
63  void CheckGraphFitOptions(Foption_t &fitOption);
64 
65 
66  void GetDrawingRange(TH1 * h1, ROOT::Fit::DataRange & range);
71 
72 
73  template <class FitObject>
74  TFitResultPtr Fit(FitObject * h1, TF1 *f1 , Foption_t & option , const ROOT::Math::MinimizerOptions & moption, const char *goption, ROOT::Fit::DataRange & range);
75 
76  template <class FitObject>
77  void StoreAndDrawFitFunction(FitObject * h1, TF1 * f1, const ROOT::Fit::DataRange & range, bool, bool, const char *goption);
78 
79  template <class FitObject>
80  double ComputeChi2(const FitObject & h1, TF1 &f1, bool useRange, bool usePL );
81 
82 
83 
84 }
85 
86 int HFit::CheckFitFunction(const TF1 * f1, int dim) {
87  // Check validity of fitted function
88  if (!f1) {
89  Error("Fit", "function may not be null pointer");
90  return -1;
91  }
92  if (f1->IsZombie()) {
93  Error("Fit", "function is zombie");
94  return -2;
95  }
96 
97  int npar = f1->GetNpar();
98  if (npar <= 0) {
99  Error("Fit", "function %s has illegal number of parameters = %d", f1->GetName(), npar);
100  return -3;
101  }
102 
103  // Check that function has same dimension as histogram
104  if (f1->GetNdim() > dim) {
105  Error("Fit","function %s dimension, %d, is greater than fit object dimension, %d",
106  f1->GetName(), f1->GetNdim(), dim);
107  return -4;
108  }
109  if (f1->GetNdim() < dim-1) {
110  Error("Fit","function %s dimension, %d, is smaller than fit object dimension -1, %d",
111  f1->GetName(), f1->GetNdim(), dim);
112  return -5;
113  }
114 
115  return 0;
116 
117 }
118 
119 
121  // get the range form the function and fill and return the DataRange object
122  Double_t fxmin, fymin, fzmin, fxmax, fymax, fzmax;
123  f1.GetRange(fxmin, fymin, fzmin, fxmax, fymax, fzmax);
124  // support only one range - so add only if was not set before
125  if (range.Size(0) == 0) range.AddRange(0,fxmin,fxmax);
126  if (range.Size(1) == 0) range.AddRange(1,fymin,fymax);
127  if (range.Size(2) == 0) range.AddRange(2,fzmin,fzmax);
128  return;
129 }
130 
131 
132 template<class FitObject>
133 TFitResultPtr HFit::Fit(FitObject * h1, TF1 *f1 , Foption_t & fitOption , const ROOT::Math::MinimizerOptions & minOption, const char *goption, ROOT::Fit::DataRange & range)
134 {
135  // perform fit of histograms, or graphs using new fitting classes
136  // use same routines for fitting both graphs and histograms
137 
138 #ifdef DEBUG
139  printf("fit function %s\n",f1->GetName() );
140 #endif
141 
142  // replacement function using new fitter
143  int hdim = HFit::GetDimension(h1);
144  int iret = HFit::CheckFitFunction(f1, hdim);
145  if (iret != 0) return iret;
146 
147 
148 
149  // integral option is not supported in this case
150  if (f1->GetNdim() < hdim ) {
151  if (fitOption.Integral) Info("Fit","Ignore Integral option. Model function dimension is less than the data object dimension");
152  if (fitOption.Like) Info("Fit","Ignore Likelihood option. Model function dimension is less than the data object dimension");
153  fitOption.Integral = 0;
154  fitOption.Like = 0;
155  }
156 
157  Int_t special = f1->GetNumber();
158  Bool_t linear = f1->IsLinear();
159  Int_t npar = f1->GetNpar();
160  if (special==299+npar) linear = kTRUE; // for polynomial functions
161  // do not use linear fitter in these case
162  if (fitOption.Bound || fitOption.Like || fitOption.Errors || fitOption.Gradient || fitOption.More || fitOption.User|| fitOption.Integral || fitOption.Minuit)
163  linear = kFALSE;
164 
165  // create an empty TFitResult
166  std::shared_ptr<TFitResult> tfr(new TFitResult() );
167  // create the fitter from an empty fit result
168  std::shared_ptr<ROOT::Fit::Fitter> fitter(new ROOT::Fit::Fitter(std::static_pointer_cast<ROOT::Fit::FitResult>(tfr) ) );
169  ROOT::Fit::FitConfig & fitConfig = fitter->Config();
170 
171  // create options
173  opt.fIntegral = fitOption.Integral;
174  opt.fUseRange = fitOption.Range;
175  opt.fExpErrors = fitOption.PChi2; // pearson chi2 with expected errors
176  if (fitOption.Like || fitOption.PChi2) opt.fUseEmpty = true; // use empty bins in log-likelihood fits
177  if (special==300) opt.fCoordErrors = false; // no need to use coordinate errors in a pol0 fit
178  if (fitOption.NoErrX) opt.fCoordErrors = false; // do not use coordinate errors when requested
179  if (fitOption.W1 ) opt.fErrors1 = true;
180  if (fitOption.W1 > 1) opt.fUseEmpty = true; // use empty bins with weight=1
181 
182  if (fitOption.BinVolume) {
183  opt.fBinVolume = true; // scale by bin volume
184  if (fitOption.BinVolume == 2) opt.fNormBinVolume = true; // scale by normalized bin volume
185  }
186 
187  if (opt.fUseRange) {
188 #ifdef DEBUG
189  printf("use range \n" );
190 #endif
191  HFit::GetFunctionRange(*f1,range);
192  }
193 #ifdef DEBUG
194  printf("range size %d\n", range.Size(0) );
195  if (range.Size(0)) {
196  double x1; double x2; range.GetRange(0,x1,x2);
197  printf(" range in x = [%f,%f] \n",x1,x2);
198  }
199 #endif
200 
201  // fill data
202  std::shared_ptr<ROOT::Fit::BinData> fitdata(new ROOT::Fit::BinData(opt,range) );
203  ROOT::Fit::FillData(*fitdata, h1, f1);
204  if (fitdata->Size() == 0 ) {
205  Warning("Fit","Fit data is empty ");
206  return -1;
207  }
208 
209 #ifdef DEBUG
210  printf("HFit:: data size is %d \n",fitdata->Size());
211  for (unsigned int i = 0; i < fitdata->Size(); ++i) {
212  if (fitdata->NDim() == 1) printf(" x[%d] = %f - value = %f \n", i,*(fitdata->Coords(i)),fitdata->Value(i) );
213  }
214 #endif
215 
216  // switch off linear fitting in case data has coordinate errors and the option is set
217  if (fitdata->GetErrorType() == ROOT::Fit::BinData::kCoordError && fitdata->Opt().fCoordErrors ) linear = false;
218  // linear fit cannot be done also in case of asymmetric errors
219  if (fitdata->GetErrorType() == ROOT::Fit::BinData::kAsymError && fitdata->Opt().fAsymErrors ) linear = false;
220 
221  // this functions use the TVirtualFitter
222  if (special != 0 && !fitOption.Bound && !linear) {
223  if (special == 100) ROOT::Fit::InitGaus (*fitdata,f1); // gaussian
224  else if (special == 110 || special == 112) ROOT::Fit::Init2DGaus(*fitdata,f1); // 2D gaussians ( xygaus or bigaus)
225  else if (special == 400) ROOT::Fit::InitGaus (*fitdata,f1); // landau (use the same)
226  else if (special == 410) ROOT::Fit::Init2DGaus(*fitdata,f1); // 2D landau (use the same)
227 
228  else if (special == 200) ROOT::Fit::InitExpo (*fitdata, f1); // exponential
229 
230  }
231 
232 
233  // set the fit function
234  // if option grad is specified use gradient
235  if ( (linear || fitOption.Gradient) )
236  fitter->SetFunction(ROOT::Math::WrappedMultiTF1(*f1));
237 #ifdef R__HAS_VECCORE
238  else if(f1->IsVectorized())
240 #endif
241  else
242  fitter->SetFunction(static_cast<const ROOT::Math::IParamMultiFunction &>(ROOT::Math::WrappedMultiTF1(*f1) ) );
243 
244  // error normalization in case of zero error in the data
245  if (fitdata->GetErrorType() == ROOT::Fit::BinData::kNoError) fitConfig.SetNormErrors(true);
246  // error normalization also in case of W or WW options (weights = 1)
247  if (fitdata->Opt().fErrors1) fitConfig.SetNormErrors(true);
248  // normalize errors also in case you are fitting a Ndim histo with a N-1 function
249  if (int(fitdata->NDim()) == hdim -1 ) fitConfig.SetNormErrors(true);
250 
251 
252  // here need to get some static extra information (like max iterations, error def, etc...)
253 
254 
255  // parameter settings and transfer the parameters values, names and limits from the functions
256  // is done automatically in the Fitter.cxx
257  for (int i = 0; i < npar; ++i) {
258  ROOT::Fit::ParameterSettings & parSettings = fitConfig.ParSettings(i);
259 
260  // check limits
261  double plow,pup;
262  f1->GetParLimits(i,plow,pup);
263  if (plow*pup != 0 && plow >= pup) { // this is a limitation - cannot fix a parameter to zero value
264  parSettings.Fix();
265  }
266  else if (plow < pup ) {
267  if (!TMath::Finite(pup) && TMath::Finite(plow) )
268  parSettings.SetLowerLimit(plow);
269  else if (!TMath::Finite(plow) && TMath::Finite(pup) )
270  parSettings.SetUpperLimit(pup);
271  else
272  parSettings.SetLimits(plow,pup);
273  }
274 
275  // set the parameter step size (by default are set to 0.3 of value)
276  // if function provides meaningful error values
277  double err = f1->GetParError(i);
278  if ( err > 0)
279  parSettings.SetStepSize(err);
280  else if (plow < pup && TMath::Finite(plow) && TMath::Finite(pup) ) { // in case of limits improve step sizes
281  double step = 0.1 * (pup - plow);
282  // check if value is not too close to limit otherwise trim value
283  if ( parSettings.Value() < pup && pup - parSettings.Value() < 2 * step )
284  step = (pup - parSettings.Value() ) / 2;
285  else if ( parSettings.Value() > plow && parSettings.Value() - plow < 2 * step )
286  step = (parSettings.Value() - plow ) / 2;
287 
288  parSettings.SetStepSize(step);
289  }
290 
291 
292  }
293 
294  // needed for setting precision ?
295  // - Compute sum of squares of errors in the bin range
296  // should maybe use stat[1] ??
297  // Double_t ey, sumw2=0;
298 // for (i=hxfirst;i<=hxlast;i++) {
299 // ey = GetBinError(i);
300 // sumw2 += ey*ey;
301 // }
302 
303 
304  // set all default minimizer options (tolerance, max iterations, etc..)
305  fitConfig.SetMinimizerOptions(minOption);
306 
307  // specific print level options
308  if (fitOption.Verbose) fitConfig.MinimizerOptions().SetPrintLevel(3);
309  if (fitOption.Quiet) fitConfig.MinimizerOptions().SetPrintLevel(0);
310 
311  // specific minimizer options depending on minimizer
312  if (linear) {
313  if (fitOption.Robust ) {
314  // robust fitting
315  std::string type = "Robust";
316  // if an h is specified print out the value adding to the type
317  if (fitOption.hRobust > 0 && fitOption.hRobust < 1.)
318  type += " (h=" + ROOT::Math::Util::ToString(fitOption.hRobust) + ")";
319  fitConfig.SetMinimizer("Linear",type.c_str());
320  fitConfig.MinimizerOptions().SetTolerance(fitOption.hRobust); // use tolerance for passing robust parameter
321  }
322  else
323  fitConfig.SetMinimizer("Linear","");
324  }
325  else {
326  if (fitOption.More) fitConfig.SetMinimizer("Minuit","MigradImproved");
327  }
328 
329 
330  // check if Error option (run Hesse and Minos) then
331  if (fitOption.Errors) {
332  // run Hesse and Minos
333  fitConfig.SetParabErrors(true);
334  fitConfig.SetMinosErrors(true);
335  }
336 
337 
338  // do fitting
339 
340 #ifdef DEBUG
341  if (fitOption.Like) printf("do likelihood fit...\n");
342  if (linear) printf("do linear fit...\n");
343  printf("do now fit...\n");
344 #endif
345 
346  bool fitok = false;
347 
348 
349  // check if can use option user
350  //typedef void (* MinuitFCN_t )(int &npar, double *gin, double &f, double *u, int flag);
351  TVirtualFitter::FCNFunc_t userFcn = 0;
352  if (fitOption.User && TVirtualFitter::GetFitter() ) {
353  userFcn = (TVirtualFitter::GetFitter())->GetFCN();
354  (TVirtualFitter::GetFitter())->SetUserFunc(f1);
355  }
356 
357 
358  if (fitOption.User && userFcn) // user provided fit objective function
359  fitok = fitter->FitFCN( userFcn );
360  else if (fitOption.Like) {// likelihood fit
361  // perform a weighted likelihood fit by applying weight correction to errors
362  bool weight = ((fitOption.Like & 2) == 2);
363  fitConfig.SetWeightCorrection(weight);
364  bool extended = ((fitOption.Like & 4 ) != 4 );
365  //if (!extended) Info("HFitImpl","Do a not -extended binned fit");
366 
367  // pass fitdata as a shared pointer so ownership is shared with Fitter and FitResult class
368  fitok = fitter->LikelihoodFit(fitdata, extended, fitOption.ExecPolicy);
369  }
370  else{ // standard least square fit
371  fitok = fitter->Fit(fitdata, fitOption.ExecPolicy);
372  }
373  if ( !fitok && !fitOption.Quiet )
374  Warning("Fit","Abnormal termination of minimization.");
375  iret |= !fitok;
376 
377 
378  const ROOT::Fit::FitResult & fitResult = fitter->Result();
379  // one could set directly the fit result in TF1
380  iret = fitResult.Status();
381  if (!fitResult.IsEmpty() ) {
382  // set in f1 the result of the fit
383  f1->SetChisquare(fitResult.Chi2() );
384  f1->SetNDF(fitResult.Ndf() );
385  f1->SetNumberFitPoints(fitdata->Size() );
386 
387  assert((Int_t)fitResult.Parameters().size() >= f1->GetNpar() );
388  f1->SetParameters( const_cast<double*>(&(fitResult.Parameters().front())));
389  if ( int( fitResult.Errors().size()) >= f1->GetNpar() )
390  f1->SetParErrors( &(fitResult.Errors().front()) );
391 
392 
393  }
394 
395 // - Store fitted function in histogram functions list and draw
396  if (!fitOption.Nostore) {
397  HFit::GetDrawingRange(h1, range);
398  HFit::StoreAndDrawFitFunction(h1, f1, range, !fitOption.Plus, !fitOption.Nograph, goption);
399  }
400 
401  // print the result
402  // if using Fitter class must be done here
403  // use old style Minuit for TMinuit and if no corrections have been applied
404  if (!fitOption.Quiet) {
405  if (fitter->GetMinimizer() && fitConfig.MinimizerType() == "Minuit" &&
406  !fitConfig.NormalizeErrors() && fitOption.Like <= 1) {
407  fitter->GetMinimizer()->PrintResults(); // use old style Minuit
408  }
409  else {
410  // print using FitResult class
411  if (fitOption.Verbose) fitResult.PrintCovMatrix(std::cout);
412  fitResult.Print(std::cout);
413  }
414  }
415 
416 
417  // store result in the backward compatible VirtualFitter
418  // in case multi-thread is not enabled
419  if (!gGlobalMutex) {
420  TVirtualFitter * lastFitter = TVirtualFitter::GetFitter();
421  TBackCompFitter * bcfitter = new TBackCompFitter(fitter, fitdata);
422  bcfitter->SetFitOption(fitOption);
423  bcfitter->SetObjectFit(h1);
424  bcfitter->SetUserFunc(f1);
426  if (userFcn) {
427  bcfitter->SetFCN(userFcn);
428  // for interpreted FCN functions
429  if (lastFitter->GetMethodCall() ) bcfitter->SetMethodCall(lastFitter->GetMethodCall() );
430  }
431 
432  // delete last fitter if it has been created here before
433  if (lastFitter) {
434  TBackCompFitter * lastBCFitter = dynamic_cast<TBackCompFitter *> (lastFitter);
435  if (lastBCFitter && lastBCFitter->TestBit(TBackCompFitter::kCanDeleteLast) )
436  delete lastBCFitter;
437  }
438  //N.B= this might create a memory leak if user does not delete the fitter they create
439  TVirtualFitter::SetFitter( bcfitter );
440  }
441 
442  // use old-style for printing the results
443  // if (fitOption.Verbose) bcfitter->PrintResults(2,0.);
444  // else if (!fitOption.Quiet) bcfitter->PrintResults(1,0.);
445 
446  if (fitOption.StoreResult)
447  {
448  TString name = "TFitResult-";
449  name = name + h1->GetName() + "-" + f1->GetName();
450  TString title = "TFitResult-";
451  title += h1->GetTitle();
452  tfr->SetName(name);
453  tfr->SetTitle(title);
454  return TFitResultPtr(tfr);
455  }
456  else
457  return TFitResultPtr(iret);
458 }
459 
460 
462  // get range from histogram and update the DataRange class
463  // if a ranges already exist in that dimension use that one
464 
465  Int_t ndim = GetDimension(h1);
466 
467  double xmin = 0, xmax = 0, ymin = 0, ymax = 0, zmin = 0, zmax = 0;
468  if (range.Size(0) == 0) {
469  TAxis & xaxis = *(h1->GetXaxis());
470  Int_t hxfirst = xaxis.GetFirst();
471  Int_t hxlast = xaxis.GetLast();
472  Double_t binwidx = xaxis.GetBinWidth(hxlast);
473  xmin = xaxis.GetBinLowEdge(hxfirst);
474  xmax = xaxis.GetBinLowEdge(hxlast) +binwidx;
475  range.AddRange(xmin,xmax);
476  }
477 
478  if (ndim > 1) {
479  if (range.Size(1) == 0) {
480  TAxis & yaxis = *(h1->GetYaxis());
481  Int_t hyfirst = yaxis.GetFirst();
482  Int_t hylast = yaxis.GetLast();
483  Double_t binwidy = yaxis.GetBinWidth(hylast);
484  ymin = yaxis.GetBinLowEdge(hyfirst);
485  ymax = yaxis.GetBinLowEdge(hylast) +binwidy;
486  range.AddRange(1,ymin,ymax);
487  }
488  }
489  if (ndim > 2) {
490  if (range.Size(2) == 0) {
491  TAxis & zaxis = *(h1->GetZaxis());
492  Int_t hzfirst = zaxis.GetFirst();
493  Int_t hzlast = zaxis.GetLast();
494  Double_t binwidz = zaxis.GetBinWidth(hzlast);
495  zmin = zaxis.GetBinLowEdge(hzfirst);
496  zmax = zaxis.GetBinLowEdge(hzlast) +binwidz;
497  range.AddRange(2,zmin,zmax);
498  }
499  }
500 #ifdef DEBUG
501  std::cout << "xmin,xmax" << xmin << " " << xmax << std::endl;
502 #endif
503 
504 }
505 
507  // get range for graph (used sub-set histogram)
508  // N.B. : this is different than in previous implementation of TGraph::Fit where range used was from xmin to xmax.
509  TH1 * h1 = gr->GetHistogram();
510  // an histogram is normally always returned for a TGraph
511  if (h1) HFit::GetDrawingRange(h1, range);
512 }
514  // get range for multi-graph (used sub-set histogram)
515  // N.B. : this is different than in previous implementation of TMultiGraph::Fit where range used was from data xmin to xmax.
516  TH1 * h1 = mg->GetHistogram();
517  if (h1) {
518  HFit::GetDrawingRange(h1, range);
519  }
520  else if (range.Size(0) == 0) {
521  // compute range from all the TGraph's belonging to the MultiGraph
522  double xmin = std::numeric_limits<double>::infinity();
523  double xmax = -std::numeric_limits<double>::infinity();
524  TIter next(mg->GetListOfGraphs() );
525  TGraph * g = 0;
526  while ( (g = (TGraph*) next() ) ) {
527  double x1 = 0, x2 = 0, y1 = 0, y2 = 0;
528  g->ComputeRange(x1,y1,x2,y2);
529  if (x1 < xmin) xmin = x1;
530  if (x2 > xmax) xmax = x2;
531  }
532  range.AddRange(xmin,xmax);
533  }
534 }
536  // get range for graph2D (used sub-set histogram)
537  // N.B. : this is different than in previous implementation of TGraph2D::Fit. There range used was always(0,0)
538  // cannot use TGraph2D::GetHistogram which makes an interpolation
539  //TH1 * h1 = gr->GetHistogram();
540  //if (h1) HFit::GetDrawingRange(h1, range);
541  // not very efficient (t.b.i.)
542  if (range.Size(0) == 0) {
543  double xmin = gr->GetXmin();
544  double xmax = gr->GetXmax();
545  range.AddRange(0,xmin,xmax);
546  }
547  if (range.Size(1) == 0) {
548  double ymin = gr->GetYmin();
549  double ymax = gr->GetYmax();
550  range.AddRange(1,ymin,ymax);
551  }
552 }
553 
555  // get range from histogram and update the DataRange class
556  // if a ranges already exist in that dimension use that one
557 
558  Int_t ndim = GetDimension(s1);
559 
560  for ( int i = 0; i < ndim; ++i ) {
561  if ( range.Size(i) == 0 ) {
562  TAxis *axis = s1->GetAxis(i);
563  range.AddRange(i, axis->GetXmin(), axis->GetXmax());
564  }
565  }
566 }
567 
568 template<class FitObject>
569 void HFit::StoreAndDrawFitFunction(FitObject * h1, TF1 * f1, const ROOT::Fit::DataRange & range, bool delOldFunction, bool drawFunction, const char *goption) {
570 // - Store fitted function in histogram functions list and draw
571 // should have separate functions for 1,2,3d ? t.b.d in case
572 
573 #ifdef DEBUG
574  std::cout <<"draw and store fit function " << f1->GetName() << std::endl;
575 #endif
576 
577 
578  Int_t ndim = GetDimension(h1);
579  double xmin = 0, xmax = 0, ymin = 0, ymax = 0, zmin = 0, zmax = 0;
580  if (range.Size(0) ) range.GetRange(0,xmin,xmax);
581  if (range.Size(1) ) range.GetRange(1,ymin,ymax);
582  if (range.Size(2) ) range.GetRange(2,zmin,zmax);
583 
584 
585 #ifdef DEBUG
586  std::cout <<"draw and store fit function " << f1->GetName()
587  << " Range in x = [ " << xmin << " , " << xmax << " ]" << std::endl;
588 #endif
589 
590  TList * funcList = h1->GetListOfFunctions();
591  if (funcList == 0){
592  Error("StoreAndDrawFitFunction","Function list has not been created - cannot store the fitted function");
593  return;
594  }
595 
596  // delete the function in the list only if
597  // the function we are fitting is not in that list
598  // If this is the case we re-use that function object and
599  // we do not create a new one (if delOldFunction is true)
600  bool reuseOldFunction = false;
601  if (delOldFunction) {
602  TIter next(funcList, kIterBackward);
603  TObject *obj;
604  while ((obj = next())) {
605  if (obj->InheritsFrom(TF1::Class())) {
606  if (obj != f1) {
607  funcList->Remove(obj);
608  delete obj;
609  }
610  else {
611  reuseOldFunction = true;
612  }
613  }
614  }
615  }
616 
617  TF1 *fnew1 = 0;
618  TF2 *fnew2 = 0;
619  TF3 *fnew3 = 0;
620 
621  // copy TF1 using TClass to avoid slicing in case of derived classes
622  if (ndim < 2) {
623  if (!reuseOldFunction) {
624  fnew1 = (TF1*)f1->IsA()->New();
625  R__ASSERT(fnew1);
626  f1->Copy(*fnew1);
627  funcList->Add(fnew1);
628  }
629  else {
630  fnew1 = f1;
631  }
632  fnew1->SetParent( h1 );
633  fnew1->SetRange(xmin,xmax);
634  fnew1->Save(xmin,xmax,0,0,0,0);
635  if (!drawFunction) fnew1->SetBit(TF1::kNotDraw);
636  fnew1->AddToGlobalList(false);
637  } else if (ndim < 3) {
638  if (!reuseOldFunction) {
639  fnew2 = (TF2*)f1->IsA()->New();
640  R__ASSERT(fnew2);
641  f1->Copy(*fnew2);
642  funcList->Add(fnew2);
643  }
644  else {
645  fnew2 = dynamic_cast<TF2*>(f1);
646  R__ASSERT(fnew2);
647  }
648  fnew2->SetRange(xmin,ymin,xmax,ymax);
649  fnew2->SetParent( h1 );
650  fnew2->Save(xmin,xmax,ymin,ymax,0,0);
651  if (!drawFunction) fnew2->SetBit(TF1::kNotDraw);
652  fnew2->AddToGlobalList(false);
653  } else {
654  if (!reuseOldFunction) {
655  fnew3 = (TF3*)f1->IsA()->New();
656  R__ASSERT(fnew3);
657  f1->Copy(*fnew3);
658  funcList->Add(fnew3);
659  }
660  else {
661  fnew2 = dynamic_cast<TF3*>(f1);
662  R__ASSERT(fnew3);
663  }
664  fnew3->SetRange(xmin,ymin,zmin,xmax,ymax,zmax);
665  fnew3->SetParent( h1 );
666  fnew3->Save(xmin,xmax,ymin,ymax,zmin,zmax);
667  if (!drawFunction) fnew3->SetBit(TF1::kNotDraw);
668  fnew3->AddToGlobalList(false);
669  }
670  if (h1->TestBit(kCanDelete)) return;
671  // draw only in case of histograms
672  if (drawFunction && ndim < 3 && h1->InheritsFrom(TH1::Class() ) ) {
673  // no need to re-draw the histogram if the histogram is already in the pad
674  // in that case the function will be just drawn (if option N is not set)
675  if (!gPad || (gPad && gPad->GetListOfPrimitives()->FindObject(h1) == NULL ) )
676  h1->Draw(goption);
677  }
678  if (gPad) gPad->Modified(); // this is not in TH1 code (needed ??)
679 
680  return;
681 }
682 
683 
684 void ROOT::Fit::FitOptionsMake(EFitObjectType type, const char *option, Foption_t &fitOption) {
685  // - Decode list of options into fitOption (used by both TGraph and TH1)
686  // works for both histograms and graph depending on the enum FitObjectType defined in HFit
689  }
690 
691  if (option == 0) return;
692  if (!option[0]) return;
693 
694  TString opt = option;
695  opt.ToUpper();
696 
697  // parse firt the specific options
698  if (type == kHistogram) {
699 
700  if (opt.Contains("WIDTH")) {
701  fitOption.BinVolume = 1; // scale content by the bin width
702  if (opt.Contains("NORMWIDTH")) {
703  // for variable bins: scale content by the bin width normalized by a reference value (typically the minimum bin)
704  // this option is for variable bin widths
705  fitOption.BinVolume = 2;
706  opt.ReplaceAll("NORMWIDTH","");
707  }
708  else
709  opt.ReplaceAll("WIDTH","");
710  }
711 
712  // if (opt.Contains("MULTIPROC")) {
713  // fitOption.ExecPolicy = ROOT::Fit::kMultiprocess;
714  // opt.ReplaceAll("MULTIPROC","");
715  // }
716 
717  if (opt.Contains("SERIAL")) {
719  opt.ReplaceAll("SERIAL","");
720  }
721 
722  if (opt.Contains("MULTITHREAD")) {
724  opt.ReplaceAll("MULTITHREAD","");
725  }
726 
727  if (opt.Contains("I")) fitOption.Integral= 1; // integral of function in the bin (no sense for graph)
728  if (opt.Contains("WW")) fitOption.W1 = 2; //all bins have weight=1, even empty bins
729  }
730 
731  // specific Graph options (need to be parsed before)
732  else if (type == kGraph) {
733  opt.ReplaceAll("ROB", "H");
734  opt.ReplaceAll("EX0", "T");
735 
736  //for robust fitting, see if # of good points is defined
737  // decode parameters for robust fitting
738  Double_t h=0;
739  if (opt.Contains("H=0.")) {
740  int start = opt.Index("H=0.");
741  int numpos = start + strlen("H=0.");
742  int numlen = 0;
743  int len = opt.Length();
744  while( (numpos+numlen<len) && isdigit(opt[numpos+numlen]) ) numlen++;
745  TString num = opt(numpos,numlen);
746  opt.Remove(start+strlen("H"),strlen("=0.")+numlen);
747  h = atof(num.Data());
748  h*=TMath::Power(10, -numlen);
749  }
750 
751  if (opt.Contains("H")) { fitOption.Robust = 1; fitOption.hRobust = h; }
752  if (opt.Contains("T")) fitOption.NoErrX = 1; // no error in X
753 
754  }
755 
756  if (opt.Contains("U")) fitOption.User = 1;
757  if (opt.Contains("Q")) fitOption.Quiet = 1;
758  if (opt.Contains("V")) {fitOption.Verbose = 1; fitOption.Quiet = 0;}
759  if (opt.Contains("L")) fitOption.Like = 1;
760  if (opt.Contains("X")) fitOption.Chi2 = 1;
761  if (opt.Contains("P")) fitOption.PChi2 = 1;
762 
763 
764  // likelihood fit options
765  if (fitOption.Like == 1) {
766  //if (opt.Contains("LL")) fitOption.Like = 2;
767  if (opt.Contains("W")){ fitOption.Like = 2; fitOption.W1=0;}// (weighted likelihood)
768  if (opt.Contains("MULTI")) {
769  if (fitOption.Like == 2) fitOption.Like = 6; // weighted multinomial
770  else fitOption.Like = 4; // multinomial likelihood fit instead of Poisson
771  opt.ReplaceAll("MULTI","");
772  }
773  // in case of histogram give precedence for likelihood options
774  if (type == kHistogram) {
775  if (fitOption.Chi2 == 1 || fitOption.PChi2 == 1)
776  Warning("Fit","Cannot use P or X option in combination of L. Ignore the chi2 option and perform a likelihood fit");
777  }
778 
779  } else {
780  if (opt.Contains("W")) fitOption.W1 = 1; // all non-empty bins have weight =1 (for chi2 fit)
781  }
782 
783  if (fitOption.PChi2 && fitOption.W1) {
784  Warning("FitOptionsMake", "Ignore option W or WW when used together with option P (Pearson chi2)");
785  fitOption.W1 = 0; // with Pearson chi2 W option is ignored
786  }
787 
788  if (opt.Contains("E")) fitOption.Errors = 1;
789  if (opt.Contains("R")) fitOption.Range = 1;
790  if (opt.Contains("G")) fitOption.Gradient= 1;
791  if (opt.Contains("M")) fitOption.More = 1;
792  if (opt.Contains("N")) fitOption.Nostore = 1;
793  if (opt.Contains("0")) fitOption.Nograph = 1;
794  if (opt.Contains("+")) fitOption.Plus = 1;
795  if (opt.Contains("B")) fitOption.Bound = 1;
796  if (opt.Contains("C")) fitOption.Nochisq = 1;
797  if (opt.Contains("F")) fitOption.Minuit = 1;
798  if (opt.Contains("S")) fitOption.StoreResult = 1;
799 
800 }
801 
803  if (foption.Like) {
804  Info("CheckGraphFitOptions","L (Log Likelihood fit) is an invalid option when fitting a graph. It is ignored");
805  foption.Like = 0;
806  }
807  if (foption.Integral) {
808  Info("CheckGraphFitOptions","I (use function integral) is an invalid option when fitting a graph. It is ignored");
809  foption.Integral = 0;
810  }
811  return;
812 }
813 
814 // implementation of unbin fit function (defined in HFitInterface)
815 
817  // do unbin fit, ownership of fitdata is passed later to the TBackFitter class
818 
819  // create a shared pointer to the fit data to managed it
820  std::shared_ptr<ROOT::Fit::UnBinData> fitdata(data);
821 
822 #ifdef DEBUG
823  printf("tree data size is %d \n",fitdata->Size());
824  for (unsigned int i = 0; i < fitdata->Size(); ++i) {
825  if (fitdata->NDim() == 1) printf(" x[%d] = %f \n", i,*(fitdata->Coords(i) ) );
826  }
827 #endif
828  if (fitdata->Size() == 0 ) {
829  Warning("Fit","Fit data is empty ");
830  return -1;
831  }
832 
833  // create an empty TFitResult
834  std::shared_ptr<TFitResult> tfr(new TFitResult() );
835  // create the fitter
836  std::shared_ptr<ROOT::Fit::Fitter> fitter(new ROOT::Fit::Fitter(tfr) );
837  ROOT::Fit::FitConfig & fitConfig = fitter->Config();
838 
839  // dimension is given by data because TF1 pointer can have wrong one
840  unsigned int dim = fitdata->NDim();
841 
842  // set the fit function
843  // if option grad is specified use gradient
844  // need to create a wrapper for an automatic normalized TF1 ???
845  if ( fitOption.Gradient ) {
846  assert ( (int) dim == fitfunc->GetNdim() );
847  fitter->SetFunction(ROOT::Math::WrappedMultiTF1(*fitfunc) );
848  }
849  else
850  fitter->SetFunction(static_cast<const ROOT::Math::IParamMultiFunction &>(ROOT::Math::WrappedMultiTF1(*fitfunc, dim) ) );
851 
852  // parameter setting is done automaticaly in the Fitter class
853  // need only to set limits
854  int npar = fitfunc->GetNpar();
855  for (int i = 0; i < npar; ++i) {
856  ROOT::Fit::ParameterSettings & parSettings = fitConfig.ParSettings(i);
857  double plow,pup;
858  fitfunc->GetParLimits(i,plow,pup);
859  // this is a limitation of TF1 interface - cannot fix a parameter to zero value
860  if (plow*pup != 0 && plow >= pup) {
861  parSettings.Fix();
862  }
863  else if (plow < pup ) {
864  if (!TMath::Finite(pup) && TMath::Finite(plow) )
865  parSettings.SetLowerLimit(plow);
866  else if (!TMath::Finite(plow) && TMath::Finite(pup) )
867  parSettings.SetUpperLimit(pup);
868  else
869  parSettings.SetLimits(plow,pup);
870  }
871 
872  // set the parameter step size (by default are set to 0.3 of value)
873  // if function provides meaningful error values
874  double err = fitfunc->GetParError(i);
875  if ( err > 0)
876  parSettings.SetStepSize(err);
877  else if (plow < pup && TMath::Finite(plow) && TMath::Finite(pup) ) { // in case of limits improve step sizes
878  double step = 0.1 * (pup - plow);
879  // check if value is not too close to limit otherwise trim value
880  if ( parSettings.Value() < pup && pup - parSettings.Value() < 2 * step )
881  step = (pup - parSettings.Value() ) / 2;
882  else if ( parSettings.Value() > plow && parSettings.Value() - plow < 2 * step )
883  step = (parSettings.Value() - plow ) / 2;
884 
885  parSettings.SetStepSize(step);
886  }
887 
888  }
889 
890  fitConfig.SetMinimizerOptions(minOption);
891 
892  if (fitOption.Verbose) fitConfig.MinimizerOptions().SetPrintLevel(3);
893  if (fitOption.Quiet) fitConfig.MinimizerOptions().SetPrintLevel(0);
894 
895  // more
896  if (fitOption.More) fitConfig.SetMinimizer("Minuit","MigradImproved");
897 
898  // chech if Minos or more options
899  if (fitOption.Errors) {
900  // run Hesse and Minos
901  fitConfig.SetParabErrors(true);
902  fitConfig.SetMinosErrors(true);
903  }
904  // use weight correction
905  if ( (fitOption.Like & 2) == 2)
906  fitConfig.SetWeightCorrection(true);
907 
908  bool extended = (fitOption.Like & 1) == 1;
909 
910  bool fitok = false;
911  fitok = fitter->LikelihoodFit(fitdata, extended, fitOption.ExecPolicy);
912  if ( !fitok && !fitOption.Quiet )
913  Warning("UnBinFit","Abnormal termination of minimization.");
914 
915  const ROOT::Fit::FitResult & fitResult = fitter->Result();
916  // one could set directly the fit result in TF1
917  int iret = fitResult.Status();
918  if (!fitResult.IsEmpty() ) {
919  // set in fitfunc the result of the fit
920  fitfunc->SetNDF(fitResult.Ndf() );
921  fitfunc->SetNumberFitPoints(fitdata->Size() );
922 
923  assert( (Int_t)fitResult.Parameters().size() >= fitfunc->GetNpar() );
924  fitfunc->SetParameters( const_cast<double*>(&(fitResult.Parameters().front())));
925  if ( int( fitResult.Errors().size()) >= fitfunc->GetNpar() )
926  fitfunc->SetParErrors( &(fitResult.Errors().front()) );
927 
928  }
929 
930  // store result in the backward compatible VirtualFitter
931  // in case not running in a multi-thread enabled mode
932  if (gGlobalMutex) {
933  TVirtualFitter * lastFitter = TVirtualFitter::GetFitter();
934  // pass ownership of Fitter and Fitdata to TBackCompFitter (fitter pointer cannot be used afterwards)
935  TBackCompFitter * bcfitter = new TBackCompFitter(fitter, fitdata);
936  // cannot use anymore now fitdata (given away ownership)
937  fitdata = 0;
938  bcfitter->SetFitOption(fitOption);
939  //bcfitter->SetObjectFit(fTree);
940  bcfitter->SetUserFunc(fitfunc);
941 
942  if (lastFitter) delete lastFitter;
943  TVirtualFitter::SetFitter( bcfitter );
944 
945  // use old-style for printing the results
946  // if (fitOption.Verbose) bcfitter->PrintResults(2,0.);
947  // else if (!fitOption.Quiet) bcfitter->PrintResults(1,0.);
948 
949  }
950  // print results
951  if (fitOption.Verbose) fitResult.PrintCovMatrix(std::cout);
952  else if (!fitOption.Quiet) fitResult.Print(std::cout);
953 
954  if (fitOption.StoreResult)
955  {
956  TString name = "TFitResult-";
957  name = name + "UnBinData-" + fitfunc->GetName();
958  TString title = "TFitResult-";
959  title += name;
960  tfr->SetName(name);
961  tfr->SetTitle(title);
962  return TFitResultPtr(tfr);
963  }
964  else
965  return TFitResultPtr(iret);
966 }
967 
968 
969 // implementations of ROOT::Fit::FitObject functions (defined in HFitInterface) in terms of the template HFit::Fit
970 
972 moption, const char *goption, ROOT::Fit::DataRange & range) {
973  // check fit options
974  // check if have weights in case of weighted likelihood
975  if ( ((foption.Like & 2) == 2) && h1->GetSumw2N() == 0) {
976  Warning("HFit::FitObject","A weighted likelihood fit is requested but histogram is not weighted - do a standard Likelihood fit");
977  foption.Like = 1;
978  }
979  // histogram fitting
980  return HFit::Fit(h1,f1,foption,moption,goption,range);
981 }
982 
983 TFitResultPtr ROOT::Fit::FitObject(TGraph * gr, TF1 *f1 , Foption_t & foption , const ROOT::Math::MinimizerOptions & moption, const char *goption, ROOT::Fit::DataRange & range) {
984  // exclude options not valid for graphs
986  // TGraph fitting
987  return HFit::Fit(gr,f1,foption,moption,goption,range);
988 }
989 
990 TFitResultPtr ROOT::Fit::FitObject(TMultiGraph * gr, TF1 *f1 , Foption_t & foption , const ROOT::Math::MinimizerOptions & moption, const char *goption, ROOT::Fit::DataRange & range) {
991  // exclude options not valid for graphs
993  // TMultiGraph fitting
994  return HFit::Fit(gr,f1,foption,moption,goption,range);
995 }
996 
997 TFitResultPtr ROOT::Fit::FitObject(TGraph2D * gr, TF1 *f1 , Foption_t & foption , const ROOT::Math::MinimizerOptions & moption, const char *goption, ROOT::Fit::DataRange & range) {
998  // exclude options not valid for graphs
1000  // TGraph2D fitting
1001  return HFit::Fit(gr,f1,foption,moption,goption,range);
1002 }
1003 
1004 TFitResultPtr ROOT::Fit::FitObject(THnBase * s1, TF1 *f1 , Foption_t & foption , const ROOT::Math::MinimizerOptions & moption, const char *goption, ROOT::Fit::DataRange & range) {
1005  // sparse histogram fitting
1006  return HFit::Fit(s1,f1,foption,moption,goption,range);
1007 }
1008 
1009 
1010 
1011 // Int_t TGraph2D::DoFit(TF2 *f2 ,Option_t *option ,Option_t *goption) {
1012 // // internal graph2D fitting methods
1013 // Foption_t fitOption;
1014 // ROOT::Fit::FitOptionsMake(option,fitOption);
1015 
1016 // // create range and minimizer options with default values
1017 // ROOT::Fit::DataRange range(2);
1018 // ROOT::Math::MinimizerOptions minOption;
1019 // return ROOT::Fit::FitObject(this, f2 , fitOption , minOption, goption, range);
1020 // }
1021 
1022 
1023 // function to compute the simple chi2 for graphs and histograms
1024 
1025 double ROOT::Fit::Chisquare(const TH1 & h1, TF1 & f1, bool useRange, bool usePL) {
1026  return HFit::ComputeChi2(h1,f1,useRange, usePL);
1027 }
1028 
1029 double ROOT::Fit::Chisquare(const TGraph & g, TF1 & f1, bool useRange) {
1030  return HFit::ComputeChi2(g,f1, useRange, false);
1031 }
1032 
1033 template<class FitObject>
1034 double HFit::ComputeChi2(const FitObject & obj, TF1 & f1, bool useRange, bool usePL ) {
1035 
1036  // implement using the fitting classes
1038  if (usePL) opt.fUseEmpty=true;
1039  ROOT::Fit::DataRange range;
1040  // get range of function
1041  if (useRange) HFit::GetFunctionRange(f1,range);
1042  // fill the data set
1043  ROOT::Fit::BinData data(opt,range);
1044  ROOT::Fit::FillData(data, &obj, &f1);
1045  if (data.Size() == 0 ) {
1046  Warning("Chisquare","data set is empty - return -1");
1047  return -1;
1048  }
1050  if (usePL) {
1051  // use the poisson log-lokelihood (Baker-Cousins chi2)
1052  ROOT::Fit::PoissonLLFunction nll(data, wf1);
1053  return 2.* nll( f1.GetParameters() ) ;
1054  }
1055  ROOT::Fit::Chi2Function chi2(data, wf1);
1056  return chi2(f1.GetParameters() );
1057 
1058 }
TGeant4Unit::mg
static constexpr double mg
Definition: TGeant4SystemOfUnits.h:210
ROOT::Fit::FitConfig::MinimizerType
const std::string & MinimizerType() const
return type of minimizer package
Definition: FitConfig.h:194
ROOT::Math::WrappedMultiTF1Templ
Class to Wrap a ROOT Function class (like TF1) in a IParamMultiFunction interface of multi-dimensions...
Definition: WrappedMultiTF1.h:57
TAxis
Class to manage histogram axis.
Definition: TAxis.h:30
ROOT::Fit::BinData::kCoordError
@ kCoordError
Definition: BinData.h:76
HFit::GetFunctionRange
void GetFunctionRange(const TF1 &f1, ROOT::Fit::DataRange &range)
Definition: HFitImpl.cxx:120
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
ROOT::Fit::FitConfig::SetWeightCorrection
void SetWeightCorrection(bool on=true)
apply the weight correction for error matric computation
Definition: FitConfig.h:239
ROOT::Fit::FitConfig::SetMinimizer
void SetMinimizer(const char *type, const char *algo=0)
set minimizer type
Definition: FitConfig.h:186
ROOT::Fit::DataOptions::fUseEmpty
bool fUseEmpty
normalize data by a normalized the bin volume (bin volume divided by a reference value)
Definition: DataOptions.h:65
TH1::GetSumw2N
virtual Int_t GetSumw2N() const
Definition: TH1.h:311
Minimizer.h
ROOT::Fit::FitConfig::SetParabErrors
void SetParabErrors(bool on=true)
set parabolic erros
Definition: FitConfig.h:233
ROOT::Fit::PoissonLikelihoodFCN
class evaluating the log likelihood for binned Poisson likelihood fits it is template to distinguish ...
Definition: PoissonLikelihoodFCN.h:58
ymax
float ymax
Definition: THbookFile.cxx:95
ROOT::Fit::FitConfig::SetMinimizerOptions
void SetMinimizerOptions(const ROOT::Math::MinimizerOptions &minopt)
set all the minimizer options using class MinimizerOptions
Definition: FitConfig.cxx:265
kTRUE
const Bool_t kTRUE
Definition: RtypesCore.h:91
HFit::GetDrawingRange
void GetDrawingRange(TH1 *h1, ROOT::Fit::DataRange &range)
Definition: HFitImpl.cxx:461
HFit::StoreAndDrawFitFunction
void StoreAndDrawFitFunction(FitObject *h1, TF1 *f1, const ROOT::Fit::DataRange &range, bool, bool, const char *goption)
Definition: HFitImpl.cxx:569
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
TObject::TestBit
R__ALWAYS_INLINE Bool_t TestBit(UInt_t f) const
Definition: TObject.h:172
ROOT::Fit::DataOptions
DataOptions : simple structure holding the options on how the data are filled.
Definition: DataOptions.h:38
TAxis::GetBinLowEdge
virtual Double_t GetBinLowEdge(Int_t bin) const
Return low edge of bin.
Definition: TAxis.cxx:518
ROOT::Fit::Chi2FCN
Chi2FCN class for binnned fits using the least square methods.
Definition: Chi2FCN.h:49
kCanDelete
@ kCanDelete
Definition: TObject.h:339
TFitResult.h
ROOT::Fit::DataOptions::fExpErrors
bool fExpErrors
use all errors equal to 1, i.e. fit without errors (default is false)
Definition: DataOptions.h:68
ROOT::Fit::ParameterSettings::SetLowerLimit
void SetLowerLimit(double low)
set a single lower limit
Definition: ParameterSettings.h:189
Warning
void Warning(const char *location, const char *msgfmt,...)
Use this function in warning situations.
Definition: TError.cxx:232
Foption_t::Errors
int Errors
Definition: Foption.h:43
ROOT::Fit::BinData
Class describing the binned data sets : vectors of x coordinates, y values and optionally error on y ...
Definition: BinData.h:62
Foption_t
Definition: Foption.h:24
ROOT::Fit::ParameterSettings::Value
double Value() const
copy constructor and assignment operators (leave them to the compiler)
Definition: ParameterSettings.h:125
ROOT::Fit::ExecutionPolicy::kMultithread
@ kMultithread
ROOT::Fit::UnBinData
Class describing the unbinned data sets (just x coordinates values) of any dimensions.
Definition: UnBinData.h:52
TF1::GetParameters
virtual Double_t * GetParameters() const
Definition: TF1.h:516
ROOT::Fit::DataRange::GetRange
void GetRange(unsigned int irange, unsigned int icoord, double &xmin, double &xmax) const
get the i-th range for given coordinate.
Definition: DataRange.h:119
TString::Data
const char * Data() const
Definition: TString.h:369
TNamed::GetTitle
virtual const char * GetTitle() const
Returns title of object.
Definition: TNamed.h:54
TGraph2D.h
TF1::GetNdim
virtual Int_t GetNdim() const
Definition: TF1.h:481
ROOT::Fit::FitOptionsMake
void FitOptionsMake(EFitObjectType type, const char *option, Foption_t &fitOption)
Decode list of options into fitOption.
Definition: HFitImpl.cxx:684
TF1::GetNpar
virtual Int_t GetNpar() const
Definition: TF1.h:477
TGraph.h
xmax
float xmax
Definition: THbookFile.cxx:95
Foption_t::Verbose
int Verbose
Definition: Foption.h:36
TBackCompFitter::SetMethodCall
virtual void SetMethodCall(TMethodCall *m)
Definition: TBackCompFitter.h:87
ROOT::Fit::FillData
void FillData(BinData &dv, const TH1 *hist, TF1 *func=0)
fill the data vector from a TH1.
Definition: HFitInterface.cxx:117
TVirtualFitter::SetFitter
static void SetFitter(TVirtualFitter *fitter, Int_t maxpar=25)
Static function to set an alternative fitter.
Definition: TVirtualFitter.cxx:257
HFitInterface.h
ROOT::Fit::Init2DGaus
void Init2DGaus(const ROOT::Fit::BinData &data, TF1 *f1)
compute initial parameter for 2D gaussian function given the fit data Set the sigma limits for zero t...
Definition: HFitInterface.cxx:395
TFitResult
Definition: TFitResult.h:34
Foption_t::Nograph
int Nograph
Definition: Foption.h:48
ROOT::Fit::DataOptions::fCoordErrors
bool fCoordErrors
use expected errors from the function and not from the data
Definition: DataOptions.h:69
PoissonLikelihoodFCN.h
TAxis::GetFirst
Int_t GetFirst() const
Return first bin on the axis i.e.
Definition: TAxis.cxx:458
TGeant4Unit::s
static constexpr double s
Definition: TGeant4SystemOfUnits.h:168
Foption_t::Nostore
int Nostore
Definition: Foption.h:47
TF2::Save
virtual void Save(Double_t xmin, Double_t xmax, Double_t ymin, Double_t ymax, Double_t zmin, Double_t zmax)
Save values of function in array fSave.
Definition: TF2.cxx:776
Foption_t::Minuit
int Minuit
Definition: Foption.h:52
HFit::GetDimension
int GetDimension(const THnBase *s1)
Definition: HFitImpl.cxx:54
ROOT::Fit::FitObject
TFitResultPtr FitObject(TH1 *h1, TF1 *f1, Foption_t &option, const ROOT::Math::MinimizerOptions &moption, const char *goption, ROOT::Fit::DataRange &range)
fitting function for a TH1 (called from TH1::Fit)
Definition: HFitImpl.cxx:971
TString::Contains
Bool_t Contains(const char *pat, ECaseCompare cmp=kExact) const
Definition: TString.h:624
TString::Length
Ssiz_t Length() const
Definition: TString.h:410
TList.h
ROOT::Fit::DataOptions::fUseRange
bool fUseRange
use empty bins (default is false) with a fixed error of 1
Definition: DataOptions.h:66
ROOT::Fit::FitConfig::SetNormErrors
void SetNormErrors(bool on=true)
set the option to normalize the error on the result according to chi2/ndf
Definition: FitConfig.h:230
TF1::SetParErrors
virtual void SetParErrors(const Double_t *errors)
Set errors for all active parameters when calling this function, the array errors must have at least ...
Definition: TF1.cxx:3505
TH1::GetListOfFunctions
TList * GetListOfFunctions() const
Definition: TH1.h:240
Foption_t::Chi2
int Chi2
Definition: Foption.h:38
TF1::kNotDraw
@ kNotDraw
Definition: TF1.h:322
FitConfig.h
TString
Definition: TString.h:136
TVirtualFitter::SetObjectFit
virtual void SetObjectFit(TObject *obj)
Definition: TVirtualFitter.h:104
ROOT::Fit::UnBinFit
TFitResultPtr UnBinFit(ROOT::Fit::UnBinData *data, TF1 *f1, Foption_t &option, const ROOT::Math::MinimizerOptions &moption)
fit an unbin data set (from tree or from histogram buffer) using a TF1 pointer and fit options.
Definition: HFitImpl.cxx:816
TF1::IsVectorized
bool IsVectorized()
Definition: TF1.h:436
gGlobalMutex
R__EXTERN TVirtualMutex * gGlobalMutex
Definition: TVirtualMutex.h:27
ROOT::Fit::FitConfig::ParSettings
const ParameterSettings & ParSettings(unsigned int i) const
get the parameter settings for the i-th parameter (const method)
Definition: FitConfig.h:81
MinimizerOptions.h
ROOT::Fit::ParameterSettings::Fix
void Fix()
fix the parameter
Definition: ParameterSettings.h:151
UnBinData.h
TVirtualFitter
Definition: TVirtualFitter.h:29
TH1::GetZaxis
TAxis * GetZaxis()
Definition: TH1.h:319
ROOT::Fit::kHistogram
@ kHistogram
Definition: HFitInterface.h:49
ROOT::Math::IParametricFunctionMultiDimTempl
IParamFunction interface (abstract class) describing multi-dimensional parameteric functions It is a ...
Definition: IParamFunction.h:104
TObject::InheritsFrom
virtual Bool_t InheritsFrom(const char *classname) const
Returns kTRUE if object inherits from class "classname".
Definition: TObject.cxx:445
Foption_t::Plus
int Plus
Definition: Foption.h:49
Foption_t::hRobust
double hRobust
Definition: Foption.h:57
h1
TH1F * h1
Definition: legend1.C:5
ROOT::Fit::FitResult::Ndf
unsigned int Ndf() const
Number of degree of freedom.
Definition: FitResult.h:168
kGraph
@ kGraph
Definition: Buttons.h:34
TBackCompFitter
Backward compatible implementation of TVirtualFitter.
Definition: TBackCompFitter.h:37
ROOT::Math::MinimizerOptions::SetPrintLevel
void SetPrintLevel(int level)
set print level
Definition: MinimizerOptions.h:131
bool
TF1::GetParLimits
virtual void GetParLimits(Int_t ipar, Double_t &parmin, Double_t &parmax) const
Return limits for parameter ipar.
Definition: TF1.cxx:1941
TString::ReplaceAll
TString & ReplaceAll(const TString &s1, const TString &s2)
Definition: TString.h:692
TF1::SetNDF
virtual void SetNDF(Int_t ndf)
Set the number of degrees of freedom ndf should be the number of points used in a fit - the number of...
Definition: TF1.cxx:3440
kIterBackward
const Bool_t kIterBackward
Definition: TCollection.h:41
x1
static const double x1[5]
Definition: RooGaussKronrodIntegrator1D.cxx:346
TMultiGraph.h
TF1::GetRange
virtual void GetRange(Double_t *xmin, Double_t *xmax) const
Return range of a generic N-D function.
Definition: TF1.cxx:2284
TF1::GetNumber
virtual Int_t GetNumber() const
Definition: TF1.h:494
Foption_t::ExecPolicy
ROOT::Fit::ExecutionPolicy ExecPolicy
Definition: Foption.h:58
TH1::GetDimension
virtual Int_t GetDimension() const
Definition: TH1.h:279
TROOT.h
TF1::SetParameters
virtual void SetParameters(const Double_t *params)
Definition: TF1.h:640
TString::ToUpper
void ToUpper()
Change string to upper case.
Definition: TString.cxx:1138
ROOT::Fit::BinData::kNoError
@ kNoError
Definition: BinData.h:76
Chi2FCN.h
ROOT::Fit::ParameterSettings::SetStepSize
void SetStepSize(double err)
set the step size
Definition: ParameterSettings.h:157
Foption_t::User
int User
Definition: Foption.h:41
TF2.h
ROOT::Fit::DataRange::AddRange
void AddRange(unsigned int icoord, double xmin, double xmax)
add a range [xmin,xmax] for the new coordinate icoord Adding a range does not delete existing one,...
Definition: DataRange.cxx:104
TGraph2D
Definition: TGraph2D.h:41
TAxis::GetXmin
Double_t GetXmin() const
Definition: TAxis.h:133
Foption_t::Gradient
int Gradient
Definition: Foption.h:46
ROOT::Fit::FitResult::PrintCovMatrix
void PrintCovMatrix(std::ostream &os) const
print error matrix and correlations
Definition: FitResult.cxx:496
Foption_t::PChi2
int PChi2
Definition: Foption.h:39
TF3.h
TFitResultPtr
Definition: TFitResultPtr.h:32
TH1::GetYaxis
TAxis * GetYaxis()
Definition: TH1.h:318
Foption_t::BinVolume
int BinVolume
Definition: Foption.h:56
xmin
float xmin
Definition: THbookFile.cxx:95
TGraph::GetHistogram
TH1F * GetHistogram() const
Returns a pointer to the histogram used to draw the axis Takes into account the two following cases.
Definition: TGraph.cxx:1482
ROOT::Fit::Fitter
Fitter class, entry point for performing all type of fits.
Definition: Fitter.h:77
h
#define h(i)
Definition: RSha256.hxx:124
TObject::SetBit
void SetBit(UInt_t f, Bool_t set)
Set or unset the user status bits as specified in f.
Definition: TObject.cxx:696
TBackCompFitter::kCanDeleteLast
@ kCanDeleteLast
Definition: TBackCompFitter.h:53
gr
TGraphErrors * gr
Definition: legend1.C:25
HFit::ComputeChi2
double ComputeChi2(const FitObject &h1, TF1 &f1, bool useRange, bool usePL)
Definition: HFitImpl.cxx:1034
ROOT::Fit::DataOptions::fIntegral
bool fIntegral
Definition: DataOptions.h:62
TString::Remove
TString & Remove(Ssiz_t pos)
Definition: TString.h:673
kFALSE
const Bool_t kFALSE
Definition: RtypesCore.h:92
s1
#define s1(x)
Definition: RSha256.hxx:109
ROOT::Fit::ParameterSettings::SetUpperLimit
void SetUpperLimit(double up)
set a single upper limit
Definition: ParameterSettings.h:182
Foption_t::Integral
int Integral
Definition: Foption.h:50
TMath::Power
LongDouble_t Power(LongDouble_t x, LongDouble_t y)
Definition: TMath.h:724
ROOT::Fit::ParameterSettings::SetLimits
void SetLimits(double low, double up)
set a double side limit, if low == up the parameter is fixed if low > up the limits are removed The c...
Definition: ParameterSettings.h:162
ROOT::Math::Util::ToString
std::string ToString(const T &val)
Utility function for conversion to strings.
Definition: Util.h:50
ROOT::Fit::FitResult
class containg the result of the fit and all the related information (fitted parameter values,...
Definition: FitResult.h:52
TF1::SetChisquare
virtual void SetChisquare(Double_t chi2)
Definition: TF1.h:608
TF1::Copy
virtual void Copy(TObject &f1) const
Copy this F1 to a new F1.
Definition: TF1.cxx:998
Foption_t::Like
int Like
Definition: Foption.h:40
ROOT::Fit::FitResult::IsEmpty
bool IsEmpty() const
True if a fit result does not exist (even invalid) with parameter values.
Definition: FitResult.h:122
ROOT::Fit::InitExpo
void InitExpo(const ROOT::Fit::BinData &data, TF1 *f1)
compute initial parameter for an exponential function given the fit data Set the constant and slope a...
Definition: HFitInterface.cxx:274
ROOT::Fit::DataRange
class describing the range in the coordinates it supports multiple range in a coordinate.
Definition: DataRange.h:45
TAxis::GetLast
Int_t GetLast() const
Return last bin on the axis i.e.
Definition: TAxis.cxx:469
ROOT::Math::MinimizerOptions
Minimizer options.
Definition: MinimizerOptions.h:42
HFit::GetDimension
int GetDimension(const TH1 *h1)
Definition: HFitImpl.cxx:50
THnBase.h
TVirtualPad.h
TF3
Definition: TF3.h:28
Foption_t::StoreResult
int StoreResult
Definition: Foption.h:55
ROOT::Fit::FitConfig::MinimizerOptions
ROOT::Math::MinimizerOptions & MinimizerOptions()
access to the minimizer control parameter (non const method)
Definition: FitConfig.h:172
TF2::SetRange
virtual void SetRange(Double_t xmin, Double_t xmax)
Initialize the upper and lower bounds to draw the function.
Definition: TF2.h:150
ROOT::Math::MinimizerOptions::SetTolerance
void SetTolerance(double tol)
set the tolerance
Definition: MinimizerOptions.h:140
HFit::CheckGraphFitOptions
void CheckGraphFitOptions(Foption_t &fitOption)
Definition: HFitImpl.cxx:802
ROOT::Fit::FitConfig
Class describing the configuration of the fit, options and parameter settings using the ROOT::Fit::Pa...
Definition: FitConfig.h:52
Foption_t::Robust
int Robust
Definition: Foption.h:54
TF1::SetNumberFitPoints
virtual void SetNumberFitPoints(Int_t npfits)
Definition: TF1.h:620
Foption_t::More
int More
Definition: Foption.h:44
ROOT::Fit::Chisquare
double Chisquare(const TH1 &h1, TF1 &f1, bool useRange, bool usePL=false)
compute the chi2 value for an histogram given a function (see TH1::Chisquare for the documentation)
Definition: HFitImpl.cxx:1025
TF1::SetRange
virtual void SetRange(Double_t xmin, Double_t xmax)
Initialize the upper and lower bounds to draw the function.
Definition: TF1.cxx:3540
ROOT::Fit::FitConfig::NormalizeErrors
bool NormalizeErrors() const
flag to check if resulting errors are be normalized according to chi2/ndf
Definition: FitConfig.h:209
ROOT::Fit::DataOptions::fErrors1
bool fErrors1
use the function range when creating the fit data (default is false)
Definition: DataOptions.h:67
ymin
float ymin
Definition: THbookFile.cxx:95
TF2
Definition: TF2.h:29
ROOT::Fit::ParameterSettings
Class, describing value, limits and step size of the parameters Provides functionality also to set/re...
Definition: ParameterSettings.h:46
TString::Index
Ssiz_t Index(const char *pat, Ssiz_t i=0, ECaseCompare cmp=kExact) const
Definition: TString.h:639
TObject::IsZombie
R__ALWAYS_INLINE Bool_t IsZombie() const
Definition: TObject.h:134
TMultiGraph
A TMultiGraph is a collection of TGraph (or derived) objects.
Definition: TMultiGraph.h:36
ROOT::Fit::DataOptions::fBinVolume
bool fBinVolume
use integral of bin content instead of bin center (default is false)
Definition: DataOptions.h:63
Foption_t::W1
int W1
Definition: Foption.h:42
ROOT::Fit::FitData::Size
unsigned int Size() const
return number of fit points
Definition: FitData.h:303
Fitter.h
ROOT::IsImplicitMTEnabled
Bool_t IsImplicitMTEnabled()
Returns true if the implicit multi-threading in ROOT is enabled.
Definition: TROOT.cxx:554
TF1::AddToGlobalList
virtual Bool_t AddToGlobalList(Bool_t on=kTRUE)
Add to global list of functions (gROOT->GetListOfFunctions() ) return previous status (true if the fu...
Definition: TF1.cxx:835
Foption_t::Bound
int Bound
Definition: Foption.h:37
f1
TF1 * f1
Definition: legend1.C:11
ROOT::Fit::DataRange::Size
unsigned int Size(unsigned int icoord=0) const
return range size for coordinate icoord (starts from zero) Size == 0 indicates no range is present [-...
Definition: DataRange.h:86
Double_t
double Double_t
Definition: RtypesCore.h:59
Foption_t::Nochisq
int Nochisq
Definition: Foption.h:51
TF1::GetParError
virtual Double_t GetParError(Int_t ipar) const
Return value of parameter number ipar.
Definition: TF1.cxx:1931
TGraph
Definition: TGraph.h:41
R__ASSERT
#define R__ASSERT(e)
Definition: TError.h:120
TF1.h
TF3::Save
virtual void Save(Double_t xmin, Double_t xmax, Double_t ymin, Double_t ymax, Double_t zmin, Double_t zmax)
Save values of function in array fSave.
Definition: TF3.cxx:530
TList::Remove
virtual TObject * Remove(TObject *obj)
Remove object from the list.
Definition: TList.cxx:821
BinData.h
TVirtualFitter::GetMethodCall
TMethodCall * GetMethodCall() const
Definition: TVirtualFitter.h:80
TFitResultPtr.h
Info
void Info(const char *location, const char *msgfmt,...)
Use this function for informational messages.
Definition: TError.cxx:221
Foption_t::Range
int Range
Definition: Foption.h:45
TBackCompFitter::SetFCN
virtual void SetFCN(void(*fcn)(Int_t &, Double_t *, Double_t &f, Double_t *, Int_t))
Override setFCN to use the Adapter to Minuit2 FCN interface To set the address of the minimization fu...
Definition: TBackCompFitter.cxx:732
TList::Add
virtual void Add(TObject *obj)
Definition: TList.h:87
TObject
Definition: TObject.h:37
TBackCompFitter.h
ROOT::Fit::FitConfig::SetMinosErrors
void SetMinosErrors(bool on=true)
set Minos erros computation to be performed after fitting
Definition: FitConfig.h:236
TH1
Definition: TH1.h:57
WrappedMultiTF1.h
ROOT::Fit::FitResult::Chi2
double Chi2() const
Chi2 fit value in case of likelihood must be computed ?
Definition: FitResult.h:165
HFit::CheckFitFunction
int CheckFitFunction(const TF1 *f1, int hdim)
Definition: HFitImpl.cxx:86
name
char name[80]
Definition: TGX11.cxx:110
HFit
Definition: HFitImpl.cxx:47
ROOT::Fit::FitResult::Status
int Status() const
minimizer status code
Definition: FitResult.h:142
x2
static const double x2[5]
Definition: RooGaussKronrodIntegrator1D.cxx:364
gPad
#define gPad
Definition: TVirtualPad.h:287
ROOT::Fit::EFitObjectType
EFitObjectType
Definition: HFitInterface.h:48
TIter
Definition: TCollection.h:233
TVirtualFitter::FCNFunc_t
void(* FCNFunc_t)(Int_t &npar, Double_t *gin, Double_t &f, Double_t *u, Int_t flag)
Definition: TVirtualFitter.h:77
TVirtualFitter::SetFitOption
virtual void SetFitOption(Foption_t option)
Definition: TVirtualFitter.h:103
Foption_t::NoErrX
int NoErrX
Definition: Foption.h:53
TNamed::GetName
virtual const char * GetName() const
Returns name of object.
Definition: TNamed.h:53
TF1::SetParent
virtual void SetParent(TObject *p=0)
Definition: TF1.h:661
type
int type
Definition: TGX11.cxx:121
TAxis::GetXmax
Double_t GetXmax() const
Definition: TAxis.h:134
ROOT::Fit::ExecutionPolicy::kSerial
@ kSerial
TVirtualFitter::GetFitter
static TVirtualFitter * GetFitter()
static: return the current Fitter
Definition: TVirtualFitter.cxx:209
TF1
1-Dim function class
Definition: TF1.h:212
Class
void Class()
Definition: Class.C:29
TH1::GetXaxis
TAxis * GetXaxis()
Get the behaviour adopted by the object about the statoverflows. See EStatOverflows for more informat...
Definition: TH1.h:317
ROOT::Fit::FitResult::Errors
const std::vector< double > & Errors() const
parameter errors (return st::vector)
Definition: FitResult.h:174
TF1::Save
virtual void Save(Double_t xmin, Double_t xmax, Double_t ymin, Double_t ymax, Double_t zmin, Double_t zmax)
Save values of function in array fSave.
Definition: TF1.cxx:3155
HFit::FitOptionsMake
void FitOptionsMake(const char *option, Foption_t &fitOption)
TH1.h
ROOT::Fit::FitResult::Parameters
const std::vector< double > & Parameters() const
parameter values (return std::vector)
Definition: FitResult.h:179
ROOT::Fit::BinData::kAsymError
@ kAsymError
Definition: BinData.h:76
ROOT::Fit::DataOptions::fNormBinVolume
bool fNormBinVolume
normalize data by the bin volume (it is used in the Poisson likelihood fits)
Definition: DataOptions.h:64
TF3::SetRange
virtual void SetRange(Double_t xmin, Double_t xmax)
Initialize the upper and lower bounds to draw the function.
Definition: TF3.h:146
TMath::Finite
Int_t Finite(Double_t x)
Check if it is finite with a mask in order to be consistent in presence of fast math.
Definition: TMath.h:760
TF1::IsLinear
virtual Bool_t IsLinear() const
Definition: TF1.h:598
TVirtualFitter::SetUserFunc
virtual void SetUserFunc(TObject *userfunc)
Definition: TVirtualFitter.h:106
THnBase
Definition: THnBase.h:43
TList
Definition: TList.h:44
TAxis::GetBinWidth
virtual Double_t GetBinWidth(Int_t bin) const
Return bin width.
Definition: TAxis.cxx:540
Foption_t::Quiet
int Quiet
Definition: Foption.h:35
TMath.h
WrappedTF1.h
int
TH1::Draw
virtual void Draw(Option_t *option="")
Draw this histogram with options.
Definition: TH1.cxx:2997
Error
void Error(const char *location, const char *msgfmt,...)
Use this function in case an error occurred.
Definition: TError.cxx:188
TError.h
ROOT::Fit::InitGaus
void InitGaus(const ROOT::Fit::BinData &data, TF1 *f1)
compute initial parameter for gaussian function given the fit data Set the sigma limits for zero top ...
Definition: HFitInterface.cxx:313
g
#define g(i)
Definition: RSha256.hxx:123