Logo ROOT   master
Reference Guide
TF1.cxx
Go to the documentation of this file.
1 // @(#)root/hist:$Id$
2 // Author: Rene Brun 18/08/95
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 "Riostream.h"
13 #include "TROOT.h"
14 #include "TMath.h"
15 #include "TF1.h"
16 #include "TH1.h"
17 #include "TGraph.h"
18 #include "TVirtualPad.h"
19 #include "TStyle.h"
20 #include "TRandom.h"
21 #include "TInterpreter.h"
22 #include "TPluginManager.h"
23 #include "TBrowser.h"
24 #include "TColor.h"
25 #include "TMethodCall.h"
26 #include "TF1Helper.h"
27 #include "TF1NormSum.h"
28 #include "TF1Convolution.h"
29 #include "TVirtualMutex.h"
30 #include "Math/WrappedFunction.h"
31 #include "Math/WrappedTF1.h"
32 #include "Math/BrentRootFinder.h"
33 #include "Math/BrentMinimizer1D.h"
34 #include "Math/BrentMethods.h"
35 #include "Math/Integrator.h"
37 #include "Math/IntegratorOptions.h"
38 #include "Math/GaussIntegrator.h"
42 #include "Math/Functor.h"
43 #include "Math/Minimizer.h"
44 #include "Math/MinimizerOptions.h"
45 #include "Math/Factory.h"
46 #include "Math/ChebyshevPol.h"
47 #include "Fit/FitResult.h"
48 // for I/O backward compatibility
49 #include "v5/TF1Data.h"
50 
51 #include "AnalyticalIntegrals.h"
52 
53 std::atomic<Bool_t> TF1::fgAbsValue(kFALSE);
55 std::atomic<Bool_t> TF1::fgAddToGlobList(kTRUE);
56 static Double_t gErrorTF1 = 0;
57 
58 using TF1Updater_t = void (*)(Int_t nobjects, TObject **from, TObject **to);
60 
61 
62 namespace {
63 struct TF1v5Convert : public TF1 {
64 public:
65  void Convert(ROOT::v5::TF1Data &from)
66  {
67  // convert old TF1 to new one
68  fNpar = from.GetNpar();
69  fNdim = from.GetNdim();
70  if (from.fType == 0) {
71  // formula functions
72  // if ndim is not 1 set xmin max to zero to avoid error in ctor
73  double xmin = from.fXmin;
74  double xmax = from.fXmax;
75  if (fNdim > 1) {
76  xmin = 0;
77  xmax = 0;
78  }
79  TF1 fnew(from.GetName(), from.GetExpFormula(), xmin, xmax);
80  if (fNdim > 1) {
81  fnew.SetRange(from.fXmin, from.fXmax);
82  }
83  fnew.Copy(*this);
84  // need to set parameter values
85  if (from.GetParameters())
86  fFormula->SetParameters(from.GetParameters());
87  } else {
88  // case of a function pointers
89  fParams = new TF1Parameters(fNpar);
90  fName = from.GetName();
91  fTitle = from.GetTitle();
92  // need to set parameter values
93  if (from.GetParameters())
94  fParams->SetParameters(from.GetParameters());
95  }
96  // copy the other data members
97  fNpx = from.fNpx;
98  fType = (EFType)from.fType;
99  fNpfits = from.fNpfits;
100  fNDF = from.fNDF;
101  fChisquare = from.fChisquare;
102  fMaximum = from.fMaximum;
103  fMinimum = from.fMinimum;
104  fXmin = from.fXmin;
105  fXmax = from.fXmax;
106 
107  if (from.fParErrors)
108  fParErrors = std::vector<Double_t>(from.fParErrors, from.fParErrors + fNpar);
109  if (from.fParMin)
110  fParMin = std::vector<Double_t>(from.fParMin, from.fParMin + fNpar);
111  if (from.fParMax)
112  fParMax = std::vector<Double_t>(from.fParMax, from.fParMax + fNpar);
113  if (from.fNsave > 0) {
114  assert(from.fSave);
115  fSave = std::vector<Double_t>(from.fSave, from.fSave + from.fNsave);
116  }
117  // set the bits
118  for (int ibit = 0; ibit < 24; ++ibit)
119  if (from.TestBit(BIT(ibit)))
120  SetBit(BIT(ibit));
121 
122  // copy the graph attributes
123  auto &fromLine = static_cast<TAttLine &>(from);
124  fromLine.Copy(*this);
125  auto &fromFill = static_cast<TAttFill &>(from);
126  fromFill.Copy(*this);
127  auto &fromMarker = static_cast<TAttMarker &>(from);
128  fromMarker.Copy(*this);
129  }
130 };
131 } // unnamed namespace
132 
133 static void R__v5TF1Updater(Int_t nobjects, TObject **from, TObject **to)
134 {
135  auto **fromv5 = (ROOT::v5::TF1Data **)from;
136  auto **target = (TF1v5Convert **)to;
137 
138  for (int i = 0; i < nobjects; ++i) {
139  if (fromv5[i] && target[i])
140  target[i]->Convert(*fromv5[i]);
141  }
142 }
143 
145 
146 ClassImp(TF1);
147 
148 // class wrapping evaluation of TF1(x) - y0
149 class GFunc {
150  const TF1 *fFunction;
151  const double fY0;
152 public:
153  GFunc(const TF1 *function , double y): fFunction(function), fY0(y) {}
154  double operator()(double x) const
155  {
156  return fFunction->Eval(x) - fY0;
157  }
158 };
159 
160 // class wrapping evaluation of -TF1(x)
161 class GInverseFunc {
162  const TF1 *fFunction;
163 public:
164  GInverseFunc(const TF1 *function): fFunction(function) {}
165 
166  double operator()(double x) const
167  {
168  return - fFunction->Eval(x);
169  }
170 };
171 // class wrapping evaluation of -TF1(x) for multi-dimension
172 class GInverseFuncNdim {
173  TF1 *fFunction;
174 public:
175  GInverseFuncNdim(TF1 *function): fFunction(function) {}
176 
177  double operator()(const double *x) const
178  {
179  return - fFunction->EvalPar(x, (Double_t *)0);
180  }
181 };
182 
183 // class wrapping function evaluation directly in 1D interface (used for integration)
184 // and implementing the methods for the momentum calculations
185 
186 class TF1_EvalWrapper : public ROOT::Math::IGenFunction {
187 public:
188  TF1_EvalWrapper(TF1 *f, const Double_t *par, bool useAbsVal, Double_t n = 1, Double_t x0 = 0) :
189  fFunc(f),
190  fPar(((par) ? par : f->GetParameters())),
191  fAbsVal(useAbsVal),
192  fN(n),
193  fX0(x0)
194  {
195  fFunc->InitArgs(fX, fPar);
196  if (par) fFunc->SetParameters(par);
197  }
198 
200  {
201  // use default copy constructor
202  TF1_EvalWrapper *f = new TF1_EvalWrapper(*this);
203  f->fFunc->InitArgs(f->fX, f->fPar);
204  return f;
205  }
206  // evaluate |f(x)|
207  Double_t DoEval(Double_t x) const
208  {
209  // use evaluation with stored parameters (i.e. pass zero)
210  fX[0] = x;
211  Double_t fval = fFunc->EvalPar(fX, 0);
212  if (fAbsVal && fval < 0) return -fval;
213  return fval;
214  }
215  // evaluate x * |f(x)|
216  Double_t EvalFirstMom(Double_t x)
217  {
218  fX[0] = x;
219  return fX[0] * TMath::Abs(fFunc->EvalPar(fX, 0));
220  }
221  // evaluate (x - x0) ^n * f(x)
222  Double_t EvalNMom(Double_t x) const
223  {
224  fX[0] = x;
225  return TMath::Power(fX[0] - fX0, fN) * TMath::Abs(fFunc->EvalPar(fX, 0));
226  }
227 
228  TF1 *fFunc;
229  mutable Double_t fX[1];
230  const double *fPar;
231  Bool_t fAbsVal;
232  Double_t fN;
233  Double_t fX0;
234 };
235 
236 ////////////////////////////////////////////////////////////////////////////////
237 /** \class TF1
238  \ingroup Hist
239  \brief 1-Dim function class
240 
241 
242 ## TF1: 1-Dim function class
243 
244 A TF1 object is a 1-Dim function defined between a lower and upper limit.
245 The function may be a simple function based on a TFormula expression or a precompiled user function.
246 The function may have associated parameters.
247 TF1 graphics function is via the TH1 and TGraph drawing functions.
248 
249 The following types of functions can be created:
250 
251 1. [Expression using variable x and no parameters]([#F1)
252 2. [Expression using variable x with parameters](#F2)
253 3. [Lambda Expression with variable x and parameters](#F3)
254 4. [A general C function with parameters](#F4)
255 5. [A general C++ function object (functor) with parameters](#F5)
256 6. [A member function with parameters of a general C++ class](#F6)
257 
258 
259 
260 ### <a name="F1"></a> 1 - Expression using variable x and no parameters
261 
262 #### Case 1: inline expression using standard C++ functions/operators
263 
264 Begin_Macro(source)
265 {
266  TF1 *fa1 = new TF1("fa1","sin(x)/x",0,10);
267  fa1->Draw();
268 }
269 End_Macro
270 
271 #### Case 2: inline expression using a ROOT function (e.g. from TMath) without parameters
272 
273 
274 Begin_Macro(source)
275 {
276  TF1 *fa2 = new TF1("fa2","TMath::DiLog(x)",0,10);
277  fa2->Draw();
278 }
279 End_Macro
280 
281 #### Case 3: inline expression using a user defined CLING function by name
282 
283 ~~~~{.cpp}
284 Double_t myFunc(double x) { return x+sin(x); }
285 ....
286 TF1 *fa3 = new TF1("fa3","myFunc(x)",-3,5);
287 fa3->Draw();
288 ~~~~
289 
290 ### <a name="F2"></a> 2 - Expression using variable x with parameters
291 
292 #### Case 1: inline expression using standard C++ functions/operators
293 
294 * Example a:
295 
296 
297 ~~~~{.cpp}
298 TF1 *fa = new TF1("fa","[0]*x*sin([1]*x)",-3,3);
299 ~~~~
300 
301 This creates a function of variable x with 2 parameters. The parameters must be initialized via:
302 
303 ~~~~{.cpp}
304  fa->SetParameter(0,value_first_parameter);
305  fa->SetParameter(1,value_second_parameter);
306 ~~~~
307 
308 
309 Parameters may be given a name:
310 
311 ~~~~{.cpp}
312  fa->SetParName(0,"Constant");
313 ~~~~
314 
315 * Example b:
316 
317 ~~~~{.cpp}
318  TF1 *fb = new TF1("fb","gaus(0)*expo(3)",0,10);
319 ~~~~
320 
321 
322 ``gaus(0)`` is a substitute for ``[0]*exp(-0.5*((x-[1])/[2])**2)`` and ``(0)`` means start numbering parameters at ``0``. ``expo(3)`` is a substitute for ``exp([3]+[4]*x)``.
323 
324 #### Case 2: inline expression using TMath functions with parameters
325 
326 ~~~~{.cpp}
327  TF1 *fb2 = new TF1("fa3","TMath::Landau(x,[0],[1],0)",-5,10);
328  fb2->SetParameters(0.2,1.3);
329  fb2->Draw();
330 ~~~~
331 
332 
333 
334 Begin_Macro
335 {
336  TCanvas *c = new TCanvas("c","c",0,0,500,300);
337  TF1 *fb2 = new TF1("fa3","TMath::Landau(x,[0],[1],0)",-5,10);
338  fb2->SetParameters(0.2,1.3); fb2->Draw();
339  return c;
340 }
341 End_Macro
342 
343 ###<a name="F3"></a> 3 - A lambda expression with variables and parameters
344 
345 \since **6.00/00:**
346 TF1 supports using lambda expressions in the formula. This allows, by using a full C++ syntax the full power of lambda
347 functions and still maintain the capability of storing the function in a file which cannot be done with
348 function pointer or lambda written not as expression, but as code (see items below).
349 
350 Example on how using lambda to define a sum of two functions.
351 Note that is necessary to provide the number of parameters
352 
353 ~~~~{.cpp}
354 TF1 f1("f1","sin(x)",0,10);
355 TF1 f2("f2","cos(x)",0,10);
356 TF1 fsum("f1","[&](double *x, double *p){ return p[0]*f1(x) + p[1]*f2(x); }",0,10,2);
357 ~~~~
358 
359 ###<a name="F4"></a> 4 - A general C function with parameters
360 
361 Consider the macro myfunc.C below:
362 
363 ~~~~{.cpp}
364  // Macro myfunc.C
365  Double_t myfunction(Double_t *x, Double_t *par)
366  {
367  Float_t xx =x[0];
368  Double_t f = TMath::Abs(par[0]*sin(par[1]*xx)/xx);
369  return f;
370  }
371  void myfunc()
372  {
373  TF1 *f1 = new TF1("myfunc",myfunction,0,10,2);
374  f1->SetParameters(2,1);
375  f1->SetParNames("constant","coefficient");
376  f1->Draw();
377  }
378  void myfit()
379  {
380  TH1F *h1=new TH1F("h1","test",100,0,10);
381  h1->FillRandom("myfunc",20000);
382  TF1 *f1 = (TF1 *)gROOT->GetFunction("myfunc");
383  f1->SetParameters(800,1);
384  h1->Fit("myfunc");
385  }
386 ~~~~
387 
388 
389 
390 In an interactive session you can do:
391 
392 ~~~~
393  Root > .L myfunc.C
394  Root > myfunc();
395  Root > myfit();
396 ~~~~
397 
398 
399 
400 TF1 objects can reference other TF1 objects of type A or B defined above. This excludes CLing or compiled functions. However, there is a restriction. A function cannot reference a basic function if the basic function is a polynomial polN.
401 
402 Example:
403 
404 ~~~~{.cpp}
405 {
406  TF1 *fcos = new TF1 ("fcos", "[0]*cos(x)", 0., 10.);
407  fcos->SetParNames( "cos");
408  fcos->SetParameter( 0, 1.1);
409 
410  TF1 *fsin = new TF1 ("fsin", "[0]*sin(x)", 0., 10.);
411  fsin->SetParNames( "sin");
412  fsin->SetParameter( 0, 2.1);
413 
414  TF1 *fsincos = new TF1 ("fsc", "fcos+fsin");
415 
416  TF1 *fs2 = new TF1 ("fs2", "fsc+fsc");
417 }
418 ~~~~
419 
420 
421 ### <a name="F5"></a> 5 - A general C++ function object (functor) with parameters
422 
423 A TF1 can be created from any C++ class implementing the operator()(double *x, double *p). The advantage of the function object is that he can have a state and reference therefore what-ever other object. In this way the user can customize his function.
424 
425 Example:
426 
427 
428 ~~~~{.cpp}
429 class MyFunctionObject {
430  public:
431  // use constructor to customize your function object
432 
433  double operator() (double *x, double *p) {
434  // function implementation using class data members
435  }
436 };
437 {
438  ....
439  MyFunctionObject fobj;
440  TF1 * f = new TF1("f",fobj,0,1,npar); // create TF1 class.
441  .....
442 }
443 ~~~~
444 
445 #### Using a lambda function as a general C++ functor object
446 
447 From C++11 we can use both std::function or even better lambda functions to create the TF1.
448 As above the lambda must have the right signature but can capture whatever we want. For example we can make
449 a TF1 from the TGraph::Eval function as shown below where we use as function parameter the graph normalization.
450 
451 ~~~~{.cpp}
452 TGraph * g = new TGraph(npointx, xvec, yvec);
453 TF1 * f = new TF1("f",[&](double*x, double *p){ return p[0]*g->Eval(x[0]); }, xmin, xmax, 1);
454 ~~~~
455 
456 
457 ### <a name="F6"></a> 6 - A member function with parameters of a general C++ class
458 
459 A TF1 can be created in this case from any member function of a class which has the signature of (double * , double *) and returning a double.
460 
461 Example:
462 
463 ~~~~{.cpp}
464 class MyFunction {
465  public:
466  ...
467  double Evaluate() (double *x, double *p) {
468  // function implementation
469  }
470 };
471 {
472  ....
473  MyFunction * fptr = new MyFunction(....); // create the user function class
474  TF1 * f = new TF1("f",fptr,&MyFunction::Evaluate,0,1,npar,"MyFunction","Evaluate"); // create TF1 class.
475 
476  .....
477 }
478 ~~~~
479 
480 See also the tutorial __math/exampleFunctor.C__ for a running example.
481 */
482 ////////////////////////////////////////////////////////////////////////////
483 
484 TF1 *TF1::fgCurrent = 0;
485 
486 
487 ////////////////////////////////////////////////////////////////////////////////
488 /// TF1 default constructor.
489 
491  TNamed(), TAttLine(), TAttFill(), TAttMarker(),
492  fXmin(0), fXmax(0), fNpar(0), fNdim(0), fType(EFType::kFormula)
493 {
494  SetFillStyle(0);
495 }
496 
497 
498 ////////////////////////////////////////////////////////////////////////////////
499 /// F1 constructor using a formula definition
500 ///
501 /// See TFormula constructor for explanation of the formula syntax.
502 ///
503 /// See tutorials: fillrandom, first, fit1, formula1, multifit
504 /// for real examples.
505 ///
506 /// Creates a function of type A or B between xmin and xmax
507 ///
508 /// if formula has the form "fffffff;xxxx;yyyy", it is assumed that
509 /// the formula string is "fffffff" and "xxxx" and "yyyy" are the
510 /// titles for the X and Y axis respectively.
511 
512 TF1::TF1(const char *name, const char *formula, Double_t xmin, Double_t xmax, EAddToList addToGlobList, bool vectorize) :
513  TNamed(name, formula), TAttLine(), TAttFill(), TAttMarker(), fType(EFType::kFormula)
514 {
515  if (xmin < xmax) {
516  fXmin = xmin;
517  fXmax = xmax;
518  } else {
519  fXmin = xmax; //when called from TF2,TF3
520  fXmax = xmin;
521  }
522  // Create rep formula (no need to add to gROOT list since we will add the TF1 object)
523 
524  // First check if we are making a convolution
525  if (TString(formula, 5) == "CONV(" && formula[strlen(formula) - 1] == ')') {
526  // Look for single ',' delimiter
527  int delimPosition = -1;
528  int parenCount = 0;
529  for (unsigned int i = 5; i < strlen(formula) - 1; i++) {
530  if (formula[i] == '(')
531  parenCount++;
532  else if (formula[i] == ')')
533  parenCount--;
534  else if (formula[i] == ',' && parenCount == 0) {
535  if (delimPosition == -1)
536  delimPosition = i;
537  else
538  Error("TF1", "CONV takes 2 arguments. Too many arguments found in : %s", formula);
539  }
540  }
541  if (delimPosition == -1)
542  Error("TF1", "CONV takes 2 arguments. Only one argument found in : %s", formula);
543 
544  // Having found the delimiter, define the first and second formulas
545  TString formula1 = TString(TString(formula)(5, delimPosition - 5));
546  TString formula2 = TString(TString(formula)(delimPosition + 1, strlen(formula) - 1 - (delimPosition + 1)));
547  // remove spaces from these formulas
548  formula1.ReplaceAll(' ', "");
549  formula2.ReplaceAll(' ', "");
550 
551  TF1 *function1 = (TF1 *)(gROOT->GetListOfFunctions()->FindObject(formula1));
552  if (function1 == nullptr)
553  function1 = new TF1((const char *)formula1, (const char *)formula1, xmin, xmax);
554  TF1 *function2 = (TF1 *)(gROOT->GetListOfFunctions()->FindObject(formula2));
555  if (function2 == nullptr)
556  function2 = new TF1((const char *)formula2, (const char *)formula2, xmin, xmax);
557 
558  // std::cout << "functions have been defined" << std::endl;
559 
560  TF1Convolution *conv = new TF1Convolution(function1, function2);
561 
562  // (note: currently ignoring `useFFT` option)
563  fNpar = conv->GetNpar();
564  fNdim = 1; // (note: may want to extend this in the future?)
565 
566  fType = EFType::kCompositionFcn;
567  fComposition = std::unique_ptr<TF1AbsComposition>(conv);
568 
569  fParams = new TF1Parameters(fNpar); // default to zeros (TF1Convolution has no GetParameters())
570  // set parameter names
571  for (int i = 0; i < fNpar; i++)
572  this->SetParName(i, conv->GetParName(i));
573  // set parameters to default values
574  int f1Npar = function1->GetNpar();
575  int f2Npar = function2->GetNpar();
576  // first, copy parameters from function1
577  for (int i = 0; i < f1Npar; i++)
578  this->SetParameter(i, function1->GetParameter(i));
579  // then, check if the "Constant" parameters were combined
580  // (this code assumes function2 has at most one parameter named "Constant")
581  if (conv->GetNpar() == f1Npar + f2Npar - 1) {
582  int cst1 = function1->GetParNumber("Constant");
583  int cst2 = function2->GetParNumber("Constant");
584  this->SetParameter(cst1, function1->GetParameter(cst1) * function2->GetParameter(cst2));
585  // and copy parameters from function2
586  for (int i = 0; i < f2Npar; i++)
587  if (i < cst2)
588  this->SetParameter(f1Npar + i, function2->GetParameter(i));
589  else if (i > cst2)
590  this->SetParameter(f1Npar + i - 1, function2->GetParameter(i));
591  } else {
592  // or if no constant, simply copy parameters from function2
593  for (int i = 0; i < f2Npar; i++)
594  this->SetParameter(i + f1Npar, function2->GetParameter(i));
595  }
596 
597  // Then check if we need NSUM syntax:
598  } else if (TString(formula, 5) == "NSUM(" && formula[strlen(formula) - 1] == ')') {
599  // using comma as delimiter
600  char delimiter = ',';
601  // first, remove "NSUM(" and ")" and spaces
602  TString formDense = TString(formula)(5,strlen(formula)-5-1);
603  formDense.ReplaceAll(' ', "");
604 
605  // make sure standard functions are defined (e.g. gaus, expo)
607 
608  // Go char-by-char to split terms and define the relevant functions
609  int parenCount = 0;
610  int termStart = 0;
611  TObjArray *newFuncs = new TObjArray();
612  newFuncs->SetOwner(kTRUE);
613  TObjArray *coeffNames = new TObjArray();
614  coeffNames->SetOwner(kTRUE);
615  TString fullFormula("");
616  for (int i = 0; i < formDense.Length(); ++i) {
617  if (formDense[i] == '(')
618  parenCount++;
619  else if (formDense[i] == ')')
620  parenCount--;
621  else if (formDense[i] == delimiter && parenCount == 0) {
622  // term goes from termStart to i
623  DefineNSUMTerm(newFuncs, coeffNames, fullFormula, formDense, termStart, i, xmin, xmax);
624  termStart = i + 1;
625  }
626  }
627  DefineNSUMTerm(newFuncs, coeffNames, fullFormula, formDense, termStart, formDense.Length(), xmin, xmax);
628 
629  TF1NormSum *normSum = new TF1NormSum(fullFormula, xmin, xmax);
630 
631  if (xmin == 0 && xmax == 1.) Info("TF1","Created TF1NormSum object using the default [0,1] range");
632 
633  fNpar = normSum->GetNpar();
634  fNdim = 1; // (note: may want to extend functionality in the future)
635 
636  fType = EFType::kCompositionFcn;
637  fComposition = std::unique_ptr<TF1AbsComposition>(normSum);
638 
639  fParams = new TF1Parameters(fNpar);
640  fParams->SetParameters(&(normSum->GetParameters())[0]); // inherit default parameters from normSum
641 
642  // Parameter names
643  for (int i = 0; i < fNpar; i++) {
644  if (coeffNames->At(i) != nullptr) {
645  TString coeffName = ((TObjString *)coeffNames->At(i))->GetString();
646  this->SetParName(i, (const char *)coeffName);
647  } else {
648  this->SetParName(i, normSum->GetParName(i));
649  }
650  }
651 
652  } else { // regular TFormula
653  fFormula = new TFormula(name, formula, false, vectorize);
654  fNpar = fFormula->GetNpar();
655  // TFormula can have dimension zero, but since this is a TF1 minimal dim is 1
656  fNdim = fFormula->GetNdim() == 0 ? 1 : fFormula->GetNdim();
657  }
658  if (fNpar) {
659  fParErrors.resize(fNpar);
660  fParMin.resize(fNpar);
661  fParMax.resize(fNpar);
662  }
663  // do we want really to have this un-documented feature where we accept cases where dim > 1
664  // by setting xmin >= xmax ??
665  if (fNdim > 1 && xmin < xmax) {
666  Error("TF1", "function: %s/%s has dimension %d instead of 1", name, formula, fNdim);
667  MakeZombie();
668  }
669 
670  DoInitialize(addToGlobList);
671 }
673  if (opt == nullptr) return TF1::EAddToList::kDefault;
674  TString option(opt);
675  option.ToUpper();
676  if (option.Contains("NL")) return TF1::EAddToList::kNo;
677  if (option.Contains("GL")) return TF1::EAddToList::kAdd;
679 }
681  if (opt == nullptr) return false;
682  TString option(opt);
683  option.ToUpper();
684  if (option.Contains("VEC")) return true;
685  return false;
686 }
687 TF1::TF1(const char *name, const char *formula, Double_t xmin, Double_t xmax, Option_t * opt) :
688 ////////////////////////////////////////////////////////////////////////////////
689 /// Same constructor as above (for TFormula based function) but passing an option strings
690 /// available options
691 /// VEC - vectorize the formula expressions (not possible for lambda based expressions)
692 /// NL - function is not stores in the global list of functions
693 /// GL - function will be always stored in the global list of functions ,
694 /// independently of the global setting of TF1::DefaultAddToGlobalList
695 ///////////////////////////////////////////////////////////////////////////////////
696  TF1(name, formula, xmin, xmax, GetGlobalListOption(opt), GetVectorizedOption(opt) )
697 {}
698 ////////////////////////////////////////////////////////////////////////////////
699 /// F1 constructor using name of an interpreted function.
700 ///
701 /// Creates a function of type C between xmin and xmax.
702 /// name is the name of an interpreted C++ function.
703 /// The function is defined with npar parameters
704 /// fcn must be a function of type:
705 ///
706 /// Double_t fcn(Double_t *x, Double_t *params)
707 ///
708 /// This constructor is called for functions of type C by the C++ interpreter.
709 ///
710 /// WARNING! A function created with this constructor cannot be Cloned.
711 
712 TF1::TF1(const char *name, Double_t xmin, Double_t xmax, Int_t npar, Int_t ndim, EAddToList addToGlobList) :
713  TF1(EFType::kInterpreted, name, xmin, xmax, npar, ndim, addToGlobList, new TF1Parameters(npar))
714 {
715  if (fName.Data()[0] == '*') { // case TF1 name starts with a *
716  Info("TF1", "TF1 has a name starting with a \'*\' - it is for saved TF1 objects in a .C file");
717  return; //case happens via SavePrimitive
718  } else if (fName.IsNull()) {
719  Error("TF1", "requires a proper function name!");
720  return;
721  }
722 
723  fMethodCall = new TMethodCall();
724  fMethodCall->InitWithPrototype(fName, "Double_t*,Double_t*");
725 
726  if (! fMethodCall->IsValid()) {
727  Error("TF1", "No function found with the signature %s(Double_t*,Double_t*)", name);
728  return;
729  }
730 }
731 
732 
733 ////////////////////////////////////////////////////////////////////////////////
734 /// Constructor using a pointer to a real function.
735 ///
736 /// \param npar is the number of free parameters used by the function
737 ///
738 /// This constructor creates a function of type C when invoked
739 /// with the normal C++ compiler.
740 ///
741 /// see test program test/stress.cxx (function stress1) for an example.
742 /// note the interface with an intermediate pointer.
743 ///
744 /// WARNING! A function created with this constructor cannot be Cloned.
745 
746 TF1::TF1(const char *name, Double_t (*fcn)(Double_t *, Double_t *), Double_t xmin, Double_t xmax, Int_t npar, Int_t ndim, EAddToList addToGlobList) :
747  TF1(EFType::kPtrScalarFreeFcn, name, xmin, xmax, npar, ndim, addToGlobList, new TF1Parameters(npar), new TF1FunctorPointerImpl<double>(ROOT::Math::ParamFunctor(fcn)))
748 {}
749 
750 ////////////////////////////////////////////////////////////////////////////////
751 /// Constructor using a pointer to real function.
752 ///
753 /// \param npar is the number of free parameters used by the function
754 ///
755 /// This constructor creates a function of type C when invoked
756 /// with the normal C++ compiler.
757 ///
758 /// see test program test/stress.cxx (function stress1) for an example.
759 /// note the interface with an intermediate pointer.
760 ///
761 /// WARNING! A function created with this constructor cannot be Cloned.
762 
763 TF1::TF1(const char *name, Double_t (*fcn)(const Double_t *, const Double_t *), Double_t xmin, Double_t xmax, Int_t npar, Int_t ndim, EAddToList addToGlobList) :
764  TF1(EFType::kPtrScalarFreeFcn, name, xmin, xmax, npar, ndim, addToGlobList, new TF1Parameters(npar), new TF1FunctorPointerImpl<double>(ROOT::Math::ParamFunctor(fcn)))
765 {}
766 
767 ////////////////////////////////////////////////////////////////////////////////
768 /// Constructor using the Functor class.
769 ///
770 /// \param xmin and
771 /// \param xmax define the plotting range of the function
772 /// \param npar is the number of free parameters used by the function
773 ///
774 /// This constructor can be used only in compiled code
775 ///
776 /// WARNING! A function created with this constructor cannot be Cloned.
777 
778 TF1::TF1(const char *name, ROOT::Math::ParamFunctor f, Double_t xmin, Double_t xmax, Int_t npar, Int_t ndim, EAddToList addToGlobList) :
779  TF1(EFType::kPtrScalarFreeFcn, name, xmin, xmax, npar, ndim, addToGlobList, new TF1Parameters(npar), new TF1FunctorPointerImpl<double>(ROOT::Math::ParamFunctor(f)))
780 {}
781 
782 ////////////////////////////////////////////////////////////////////////////////
783 /// Common initialization of the TF1. Add to the global list and
784 /// set the default style
785 
786 void TF1::DoInitialize(EAddToList addToGlobalList)
787 {
788  // add to global list of functions if default adding is on OR if bit is set
789  bool doAdd = ((addToGlobalList == EAddToList::kDefault && fgAddToGlobList)
790  || addToGlobalList == EAddToList::kAdd);
791  if (doAdd && gROOT) {
794  // Store formula in linked list of formula in ROOT
795  TF1 *f1old = (TF1 *)gROOT->GetListOfFunctions()->FindObject(fName);
796  if (f1old) {
797  gROOT->GetListOfFunctions()->Remove(f1old);
798  // We removed f1old from the list, it is not longer global.
799  // (See TF1::AddToGlobalList which requires this flag to be correct).
800  f1old->SetBit(kNotGlobal, kTRUE);
801  }
802  gROOT->GetListOfFunctions()->Add(this);
803  } else
805 
806  if (gStyle) {
810  }
811  SetFillStyle(0);
812 }
813 
814 ////////////////////////////////////////////////////////////////////////////////
815 /// Static method to add/avoid to add automatically functions to the global list (gROOT->GetListOfFunctions() )
816 /// After having called this static method, all the functions created afterwards will follow the
817 /// desired behaviour.
818 ///
819 /// By default the functions are added automatically
820 /// It returns the previous status (true if the functions are added automatically)
821 
823 {
824  return fgAddToGlobList.exchange(on);
825 }
826 
827 ////////////////////////////////////////////////////////////////////////////////
828 /// Add to global list of functions (gROOT->GetListOfFunctions() )
829 /// return previous status (true if the function was already in the list false if not)
830 
832 {
833  if (!gROOT) return false;
834 
835  bool prevStatus = !TestBit(kNotGlobal);
836  if (on) {
837  if (prevStatus) {
839  assert(gROOT->GetListOfFunctions()->FindObject(this) != nullptr);
840  return on; // do nothing
841  }
842  // do I need to delete previous one with the same name ???
843  //TF1 * old = dynamic_cast<TF1*>( gROOT->GetListOfFunctions()->FindObject(GetName()) );
844  //if (old) { gROOT->GetListOfFunctions()->Remove(old); old->SetBit(kNotGlobal, kTRUE); }
846  gROOT->GetListOfFunctions()->Add(this);
848  } else if (prevStatus) {
849  // if previous status was on and now is off we need to remove the function
852  TF1 *old = dynamic_cast<TF1 *>(gROOT->GetListOfFunctions()->FindObject(GetName()));
853  if (!old) {
854  Warning("AddToGlobalList", "Function is supposed to be in the global list but it is not present");
855  return kFALSE;
856  }
857  gROOT->GetListOfFunctions()->Remove(this);
858  }
859  return prevStatus;
860 }
861 
862 ////////////////////////////////////////////////////////////////////////////////
863 /// Helper functions for NSUM parsing
864 
865 // Defines the formula that a given term uses, if not already defined,
866 // and appends "sanitized" formula to `fullFormula` string
867 void TF1::DefineNSUMTerm(TObjArray *newFuncs, TObjArray *coeffNames, TString &fullFormula, TString &formula,
868  int termStart, int termEnd, Double_t xmin, Double_t xmax)
869 {
870  TString originalTerm = formula(termStart, termEnd-termStart);
871  int coeffLength = TermCoeffLength(originalTerm);
872  if (coeffLength != -1)
873  termStart += coeffLength + 1;
874 
875  // `originalFunc` is the real formula and `cleanedFunc` is the
876  // sanitized version that will not confuse the TF1NormSum
877  // constructor
878  TString originalFunc = formula(termStart, termEnd-termStart);
879  TString cleanedFunc = TString(formula(termStart, termEnd-termStart))
880  .ReplaceAll('+', "<plus>")
881  .ReplaceAll('*',"<times>");
882 
883  // define function (if necessary)
884  if (!gROOT->GetListOfFunctions()->FindObject(cleanedFunc))
885  newFuncs->Add(new TF1(cleanedFunc, originalFunc, xmin, xmax));
886 
887  // append sanitized term to `fullFormula`
888  if (fullFormula.Length() != 0)
889  fullFormula.Append('+');
890 
891  // include numerical coefficient
892  if (coeffLength != -1 && originalTerm[0] != '[')
893  fullFormula.Append(originalTerm(0, coeffLength+1));
894 
895  // add coefficient name
896  if (coeffLength != -1 && originalTerm[0] == '[')
897  coeffNames->Add(new TObjString(TString(originalTerm(1,coeffLength-2))));
898  else
899  coeffNames->Add(nullptr);
900 
901  fullFormula.Append(cleanedFunc);
902 }
903 
904 
905 // Returns length of coeff at beginning of a given term, not counting the '*'
906 // Returns -1 if no coeff found
907 // Coeff can be either a number or parameter name
909  int firstAsterisk = term.First('*');
910  if (firstAsterisk == -1) // no asterisk found
911  return -1;
912 
913  if (TString(term(0,firstAsterisk)).IsFloat())
914  return firstAsterisk;
915 
916  if (term[0] == '[' && term[firstAsterisk-1] == ']'
917  && TString(term(1,firstAsterisk-2)).IsAlnum())
918  return firstAsterisk;
919 
920  return -1;
921 }
922 
923 ////////////////////////////////////////////////////////////////////////////////
924 /// Operator =
925 
926 TF1 &TF1::operator=(const TF1 &rhs)
927 {
928  if (this != &rhs) {
929  rhs.Copy(*this);
930  }
931  return *this;
932 }
933 
934 
935 ////////////////////////////////////////////////////////////////////////////////
936 /// TF1 default destructor.
937 
939 {
940  if (fHistogram) delete fHistogram;
941  if (fMethodCall) delete fMethodCall;
942 
943  // this was before in TFormula destructor
944  {
946  if (gROOT) gROOT->GetListOfFunctions()->Remove(this);
947  }
948 
949  if (fParent) fParent->RecursiveRemove(this);
950 
951  if (fFormula) delete fFormula;
952  if (fParams) delete fParams;
953  if (fFunctor) delete fFunctor;
954 }
955 
956 
957 ////////////////////////////////////////////////////////////////////////////////
958 
959 TF1::TF1(const TF1 &f1) :
961  fXmin(0), fXmax(0), fNpar(0), fNdim(0), fType(EFType::kFormula)
962 {
963  ((TF1 &)f1).Copy(*this);
964 }
965 
966 
967 ////////////////////////////////////////////////////////////////////////////////
968 /// Static function: set the fgAbsValue flag.
969 /// By default TF1::Integral uses the original function value to compute the integral
970 /// However, TF1::Moment, CentralMoment require to compute the integral
971 /// using the absolute value of the function.
972 
974 {
975  fgAbsValue = flag;
976 }
977 
978 
979 ////////////////////////////////////////////////////////////////////////////////
980 /// Browse.
981 
983 {
984  Draw(b ? b->GetDrawOption() : "");
985  gPad->Update();
986 }
987 
988 
989 ////////////////////////////////////////////////////////////////////////////////
990 /// Copy this F1 to a new F1.
991 /// Note that the cached integral with its related arrays are not copied
992 /// (they are also set as transient data members)
993 
994 void TF1::Copy(TObject &obj) const
995 {
996  delete((TF1 &)obj).fHistogram;
997  delete((TF1 &)obj).fMethodCall;
998 
999  TNamed::Copy((TF1 &)obj);
1000  TAttLine::Copy((TF1 &)obj);
1001  TAttFill::Copy((TF1 &)obj);
1002  TAttMarker::Copy((TF1 &)obj);
1003  ((TF1 &)obj).fXmin = fXmin;
1004  ((TF1 &)obj).fXmax = fXmax;
1005  ((TF1 &)obj).fNpx = fNpx;
1006  ((TF1 &)obj).fNpar = fNpar;
1007  ((TF1 &)obj).fNdim = fNdim;
1008  ((TF1 &)obj).fType = fType;
1009  ((TF1 &)obj).fChisquare = fChisquare;
1010  ((TF1 &)obj).fNpfits = fNpfits;
1011  ((TF1 &)obj).fNDF = fNDF;
1012  ((TF1 &)obj).fMinimum = fMinimum;
1013  ((TF1 &)obj).fMaximum = fMaximum;
1014 
1015  ((TF1 &)obj).fParErrors = fParErrors;
1016  ((TF1 &)obj).fParMin = fParMin;
1017  ((TF1 &)obj).fParMax = fParMax;
1018  ((TF1 &)obj).fParent = fParent;
1019  ((TF1 &)obj).fSave = fSave;
1020  ((TF1 &)obj).fHistogram = 0;
1021  ((TF1 &)obj).fMethodCall = 0;
1022  ((TF1 &)obj).fNormalized = fNormalized;
1023  ((TF1 &)obj).fNormIntegral = fNormIntegral;
1024  ((TF1 &)obj).fFormula = 0;
1025 
1026  if (fFormula) assert(fFormula->GetNpar() == fNpar);
1027 
1028  if (fMethodCall) {
1029  // use copy-constructor of TMethodCall
1030  if (((TF1 &)obj).fMethodCall) delete((TF1 &)obj).fMethodCall;
1032 // m->InitWithPrototype(fMethodCall->GetMethodName(),fMethodCall->GetProto());
1033  ((TF1 &)obj).fMethodCall = m;
1034  }
1035  if (fFormula) {
1036  TFormula *formulaToCopy = ((TF1 &)obj).fFormula;
1037  if (formulaToCopy) delete formulaToCopy;
1038  formulaToCopy = new TFormula();
1039  fFormula->Copy(*formulaToCopy);
1040  ((TF1 &)obj).fFormula = formulaToCopy;
1041  }
1042  if (fParams) {
1043  TF1Parameters *paramsToCopy = ((TF1 &)obj).fParams;
1044  if (paramsToCopy) *paramsToCopy = *fParams;
1045  else ((TF1 &)obj).fParams = new TF1Parameters(*fParams);
1046  }
1047  if (fFunctor) {
1048  // use clone of TF1FunctorPointer
1049  if (((TF1 &)obj).fFunctor) delete((TF1 &)obj).fFunctor;
1050  ((TF1 &)obj).fFunctor = fFunctor->Clone();
1051  }
1052 
1053  if (fComposition) {
1054  TF1AbsComposition *comp = (TF1AbsComposition *)fComposition->IsA()->New();
1055  fComposition->Copy(*comp);
1056  ((TF1 &)obj).fComposition = std::unique_ptr<TF1AbsComposition>(comp);
1057  }
1058 }
1059 
1060 
1061 ////////////////////////////////////////////////////////////////////////////////
1062 /// Make a complete copy of the underlying object. If 'newname' is set,
1063 /// the copy's name will be set to that name.
1064 
1065 TObject* TF1::Clone(const char* newname) const
1066 {
1067 
1068  TF1* obj = (TF1*) TNamed::Clone(newname);
1069 
1070  if (fHistogram) {
1071  obj->fHistogram = (TH1*)fHistogram->Clone();
1072  obj->fHistogram->SetDirectory(0);
1073  }
1074 
1075  return obj;
1076 }
1077 
1078 
1079 ////////////////////////////////////////////////////////////////////////////////
1080 /// Returns the first derivative of the function at point x,
1081 /// computed by Richardson's extrapolation method (use 2 derivative estimates
1082 /// to compute a third, more accurate estimation)
1083 /// first, derivatives with steps h and h/2 are computed by central difference formulas
1084 /// \f[
1085 /// D(h) = \frac{f(x+h) - f(x-h)}{2h}
1086 /// \f]
1087 /// the final estimate
1088 /// \f[
1089 /// D = \frac{4D(h/2) - D(h)}{3}
1090 /// \f]
1091 /// "Numerical Methods for Scientists and Engineers", H.M.Antia, 2nd edition"
1092 ///
1093 /// if the argument params is null, the current function parameters are used,
1094 /// otherwise the parameters in params are used.
1095 ///
1096 /// the argument eps may be specified to control the step size (precision).
1097 /// the step size is taken as eps*(xmax-xmin).
1098 /// the default value (0.001) should be good enough for the vast majority
1099 /// of functions. Give a smaller value if your function has many changes
1100 /// of the second derivative in the function range.
1101 ///
1102 /// Getting the error via TF1::DerivativeError:
1103 /// (total error = roundoff error + interpolation error)
1104 /// the estimate of the roundoff error is taken as follows:
1105 /// \f[
1106 /// err = k\sqrt{f(x)^{2} + x^{2}deriv^{2}}\sqrt{\sum ai^{2}},
1107 /// \f]
1108 /// where k is the double precision, ai are coefficients used in
1109 /// central difference formulas
1110 /// interpolation error is decreased by making the step size h smaller.
1111 ///
1112 /// \author Anna Kreshuk
1113 
1115 {
1116  if (GetNdim() > 1) {
1117  Warning("Derivative", "Function dimension is larger than one");
1118  }
1119 
1121  double xmin, xmax;
1122  GetRange(xmin, xmax);
1123  // this is not optimal (should be used the average x instead of the range)
1124  double h = eps * std::abs(xmax - xmin);
1125  if (h <= 0) h = 0.001;
1126  double der = 0;
1127  if (params) {
1128  ROOT::Math::WrappedTF1 wtf(*(const_cast<TF1 *>(this)));
1129  wtf.SetParameters(params);
1130  der = rd.Derivative1(wtf, x, h);
1131  } else {
1132  // no need to set parameters used a non-parametric wrapper to avoid allocating
1133  // an array with parameter values
1135  der = rd.Derivative1(wf, x, h);
1136  }
1137 
1138  gErrorTF1 = rd.Error();
1139  return der;
1140 
1141 }
1142 
1143 
1144 ////////////////////////////////////////////////////////////////////////////////
1145 /// Returns the second derivative of the function at point x,
1146 /// computed by Richardson's extrapolation method (use 2 derivative estimates
1147 /// to compute a third, more accurate estimation)
1148 /// first, derivatives with steps h and h/2 are computed by central difference formulas
1149 /// \f[
1150 /// D(h) = \frac{f(x+h) - 2f(x) + f(x-h)}{h^{2}}
1151 /// \f]
1152 /// the final estimate
1153 /// \f[
1154 /// D = \frac{4D(h/2) - D(h)}{3}
1155 /// \f]
1156 /// "Numerical Methods for Scientists and Engineers", H.M.Antia, 2nd edition"
1157 ///
1158 /// if the argument params is null, the current function parameters are used,
1159 /// otherwise the parameters in params are used.
1160 ///
1161 /// the argument eps may be specified to control the step size (precision).
1162 /// the step size is taken as eps*(xmax-xmin).
1163 /// the default value (0.001) should be good enough for the vast majority
1164 /// of functions. Give a smaller value if your function has many changes
1165 /// of the second derivative in the function range.
1166 ///
1167 /// Getting the error via TF1::DerivativeError:
1168 /// (total error = roundoff error + interpolation error)
1169 /// the estimate of the roundoff error is taken as follows:
1170 /// \f[
1171 /// err = k\sqrt{f(x)^{2} + x^{2}deriv^{2}}\sqrt{\sum ai^{2}},
1172 /// \f]
1173 /// where k is the double precision, ai are coefficients used in
1174 /// central difference formulas
1175 /// interpolation error is decreased by making the step size h smaller.
1176 ///
1177 /// \author Anna Kreshuk
1178 
1180 {
1181  if (GetNdim() > 1) {
1182  Warning("Derivative2", "Function dimension is larger than one");
1183  }
1184 
1186  double xmin, xmax;
1187  GetRange(xmin, xmax);
1188  // this is not optimal (should be used the average x instead of the range)
1189  double h = eps * std::abs(xmax - xmin);
1190  if (h <= 0) h = 0.001;
1191  double der = 0;
1192  if (params) {
1193  ROOT::Math::WrappedTF1 wtf(*(const_cast<TF1 *>(this)));
1194  wtf.SetParameters(params);
1195  der = rd.Derivative2(wtf, x, h);
1196  } else {
1197  // no need to set parameters used a non-parametric wrapper to avoid allocating
1198  // an array with parameter values
1200  der = rd.Derivative2(wf, x, h);
1201  }
1202 
1203  gErrorTF1 = rd.Error();
1204 
1205  return der;
1206 }
1207 
1208 
1209 ////////////////////////////////////////////////////////////////////////////////
1210 /// Returns the third derivative of the function at point x,
1211 /// computed by Richardson's extrapolation method (use 2 derivative estimates
1212 /// to compute a third, more accurate estimation)
1213 /// first, derivatives with steps h and h/2 are computed by central difference formulas
1214 /// \f[
1215 /// D(h) = \frac{f(x+2h) - 2f(x+h) + 2f(x-h) - f(x-2h)}{2h^{3}}
1216 /// \f]
1217 /// the final estimate
1218 /// \f[
1219 /// D = \frac{4D(h/2) - D(h)}{3}
1220 /// \f]
1221 /// "Numerical Methods for Scientists and Engineers", H.M.Antia, 2nd edition"
1222 ///
1223 /// if the argument params is null, the current function parameters are used,
1224 /// otherwise the parameters in params are used.
1225 ///
1226 /// the argument eps may be specified to control the step size (precision).
1227 /// the step size is taken as eps*(xmax-xmin).
1228 /// the default value (0.001) should be good enough for the vast majority
1229 /// of functions. Give a smaller value if your function has many changes
1230 /// of the second derivative in the function range.
1231 ///
1232 /// Getting the error via TF1::DerivativeError:
1233 /// (total error = roundoff error + interpolation error)
1234 /// the estimate of the roundoff error is taken as follows:
1235 /// \f[
1236 /// err = k\sqrt{f(x)^{2} + x^{2}deriv^{2}}\sqrt{\sum ai^{2}},
1237 /// \f]
1238 /// where k is the double precision, ai are coefficients used in
1239 /// central difference formulas
1240 /// interpolation error is decreased by making the step size h smaller.
1241 ///
1242 /// \author Anna Kreshuk
1243 
1245 {
1246  if (GetNdim() > 1) {
1247  Warning("Derivative3", "Function dimension is larger than one");
1248  }
1249 
1251  double xmin, xmax;
1252  GetRange(xmin, xmax);
1253  // this is not optimal (should be used the average x instead of the range)
1254  double h = eps * std::abs(xmax - xmin);
1255  if (h <= 0) h = 0.001;
1256  double der = 0;
1257  if (params) {
1258  ROOT::Math::WrappedTF1 wtf(*(const_cast<TF1 *>(this)));
1259  wtf.SetParameters(params);
1260  der = rd.Derivative3(wtf, x, h);
1261  } else {
1262  // no need to set parameters used a non-parametric wrapper to avoid allocating
1263  // an array with parameter values
1265  der = rd.Derivative3(wf, x, h);
1266  }
1267 
1268  gErrorTF1 = rd.Error();
1269  return der;
1270 
1271 }
1272 
1273 
1274 ////////////////////////////////////////////////////////////////////////////////
1275 /// Static function returning the error of the last call to the of Derivative's
1276 /// functions
1277 
1279 {
1280  return gErrorTF1;
1281 }
1282 
1283 
1284 ////////////////////////////////////////////////////////////////////////////////
1285 /// Compute distance from point px,py to a function.
1286 ///
1287 /// Compute the closest distance of approach from point px,py to this
1288 /// function. The distance is computed in pixels units.
1289 ///
1290 /// Note that px is called with a negative value when the TF1 is in
1291 /// TGraph or TH1 list of functions. In this case there is no point
1292 /// looking at the histogram axis.
1293 
1295 {
1296  if (!fHistogram) return 9999;
1297  Int_t distance = 9999;
1298  if (px >= 0) {
1299  distance = fHistogram->DistancetoPrimitive(px, py);
1300  if (distance <= 1) return distance;
1301  } else {
1302  px = -px;
1303  }
1304 
1305  Double_t xx[1];
1306  Double_t x = gPad->AbsPixeltoX(px);
1307  xx[0] = gPad->PadtoX(x);
1308  if (xx[0] < fXmin || xx[0] > fXmax) return distance;
1309  Double_t fval = Eval(xx[0]);
1310  Double_t y = gPad->YtoPad(fval);
1311  Int_t pybin = gPad->YtoAbsPixel(y);
1312  return TMath::Abs(py - pybin);
1313 }
1314 
1315 
1316 ////////////////////////////////////////////////////////////////////////////////
1317 /// Draw this function with its current attributes.
1318 ///
1319 /// Possible option values are:
1320 ///
1321 /// option | description
1322 /// -------|----------------------------------------
1323 /// "SAME" | superimpose on top of existing picture
1324 /// "L" | connect all computed points with a straight line
1325 /// "C" | connect all computed points with a smooth curve
1326 /// "FC" | draw a fill area below a smooth curve
1327 ///
1328 /// Note that the default value is "L". Therefore to draw on top
1329 /// of an existing picture, specify option "LSAME"
1330 ///
1331 /// NB. You must use DrawCopy if you want to draw several times the same
1332 /// function in the current canvas.
1333 
1334 void TF1::Draw(Option_t *option)
1335 {
1336  TString opt = option;
1337  opt.ToLower();
1338  if (gPad && !opt.Contains("same")) gPad->Clear();
1339 
1340  AppendPad(option);
1341 
1342  gPad->IncrementPaletteColor(1, opt);
1343 }
1344 
1345 
1346 ////////////////////////////////////////////////////////////////////////////////
1347 /// Draw a copy of this function with its current attributes.
1348 ///
1349 /// This function MUST be used instead of Draw when you want to draw
1350 /// the same function with different parameters settings in the same canvas.
1351 ///
1352 /// Possible option values are:
1353 ///
1354 /// option | description
1355 /// -------|----------------------------------------
1356 /// "SAME" | superimpose on top of existing picture
1357 /// "L" | connect all computed points with a straight line
1358 /// "C" | connect all computed points with a smooth curve
1359 /// "FC" | draw a fill area below a smooth curve
1360 ///
1361 /// Note that the default value is "L". Therefore to draw on top
1362 /// of an existing picture, specify option "LSAME"
1363 
1364 TF1 *TF1::DrawCopy(Option_t *option) const
1365 {
1366  TF1 *newf1 = (TF1 *)this->IsA()->New();
1367  Copy(*newf1);
1368  newf1->AppendPad(option);
1369  newf1->SetBit(kCanDelete);
1370  return newf1;
1371 }
1372 
1373 
1374 ////////////////////////////////////////////////////////////////////////////////
1375 /// Draw derivative of this function
1376 ///
1377 /// An intermediate TGraph object is built and drawn with option.
1378 /// The function returns a pointer to the TGraph object. Do:
1379 ///
1380 /// TGraph *g = (TGraph*)myfunc.DrawDerivative(option);
1381 ///
1382 /// The resulting graph will be drawn into the current pad.
1383 /// If this function is used via the context menu, it recommended
1384 /// to create a new canvas/pad before invoking this function.
1385 
1387 {
1388  TVirtualPad *pad = gROOT->GetSelectedPad();
1389  TVirtualPad *padsav = gPad;
1390  if (pad) pad->cd();
1391 
1392  TGraph *gr = new TGraph(this, "d");
1393  gr->Draw(option);
1394  if (padsav) padsav->cd();
1395  return gr;
1396 }
1397 
1398 
1399 ////////////////////////////////////////////////////////////////////////////////
1400 /// Draw integral of this function
1401 ///
1402 /// An intermediate TGraph object is built and drawn with option.
1403 /// The function returns a pointer to the TGraph object. Do:
1404 ///
1405 /// TGraph *g = (TGraph*)myfunc.DrawIntegral(option);
1406 ///
1407 /// The resulting graph will be drawn into the current pad.
1408 /// If this function is used via the context menu, it recommended
1409 /// to create a new canvas/pad before invoking this function.
1410 
1412 {
1413  TVirtualPad *pad = gROOT->GetSelectedPad();
1414  TVirtualPad *padsav = gPad;
1415  if (pad) pad->cd();
1416 
1417  TGraph *gr = new TGraph(this, "i");
1418  gr->Draw(option);
1419  if (padsav) padsav->cd();
1420  return gr;
1421 }
1422 
1423 
1424 ////////////////////////////////////////////////////////////////////////////////
1425 /// Draw function between xmin and xmax.
1426 
1428 {
1429 // //if(Compile(formula)) return ;
1430  SetRange(xmin, xmax);
1431 
1432  Draw(option);
1433 }
1434 
1435 
1436 ////////////////////////////////////////////////////////////////////////////////
1437 /// Evaluate this function.
1438 ///
1439 /// Computes the value of this function (general case for a 3-d function)
1440 /// at point x,y,z.
1441 /// For a 1-d function give y=0 and z=0
1442 /// The current value of variables x,y,z is passed through x, y and z.
1443 /// The parameters used will be the ones in the array params if params is given
1444 /// otherwise parameters will be taken from the stored data members fParams
1445 
1447 {
1448  if (fType == EFType::kFormula) return fFormula->Eval(x, y, z, t);
1449 
1450  Double_t xx[4] = {x, y, z, t};
1451  Double_t *pp = (Double_t *)fParams->GetParameters();
1452  // if (fType == EFType::kInterpreted)((TF1 *)this)->InitArgs(xx, pp);
1453  return ((TF1 *)this)->EvalPar(xx, pp);
1454 }
1455 
1456 
1457 ////////////////////////////////////////////////////////////////////////////////
1458 /// Evaluate function with given coordinates and parameters.
1459 ///
1460 /// Compute the value of this function at point defined by array x
1461 /// and current values of parameters in array params.
1462 /// If argument params is omitted or equal 0, the internal values
1463 /// of parameters (array fParams) will be used instead.
1464 /// For a 1-D function only x[0] must be given.
1465 /// In case of a multi-dimensional function, the arrays x must be
1466 /// filled with the corresponding number of dimensions.
1467 ///
1468 /// WARNING. In case of an interpreted function (fType=2), it is the
1469 /// user's responsibility to initialize the parameters via InitArgs
1470 /// before calling this function.
1471 /// InitArgs should be called at least once to specify the addresses
1472 /// of the arguments x and params.
1473 /// InitArgs should be called every time these addresses change.
1474 
1475 Double_t TF1::EvalPar(const Double_t *x, const Double_t *params)
1476 {
1477  //fgCurrent = this;
1478 
1479  if (fType == EFType::kFormula) {
1480  assert(fFormula);
1481 
1482  if (fNormalized && fNormIntegral != 0)
1483  return fFormula->EvalPar(x, params) / fNormIntegral;
1484  else
1485  return fFormula->EvalPar(x, params);
1486  }
1487  Double_t result = 0;
1488  if (fType == EFType::kPtrScalarFreeFcn || fType == EFType::kTemplScalar) {
1489  if (fFunctor) {
1490  assert(fParams);
1491  if (params) result = ((TF1FunctorPointerImpl<Double_t> *)fFunctor)->fImpl((Double_t *)x, (Double_t *)params);
1492  else result = ((TF1FunctorPointerImpl<Double_t> *)fFunctor)->fImpl((Double_t *)x, (Double_t *)fParams->GetParameters());
1493 
1494  } else result = GetSave(x);
1495 
1496  if (fNormalized && fNormIntegral != 0)
1497  result = result / fNormIntegral;
1498 
1499  return result;
1500  }
1501  if (fType == EFType::kInterpreted) {
1502  if (fMethodCall) fMethodCall->Execute(result);
1503  else result = GetSave(x);
1504 
1505  if (fNormalized && fNormIntegral != 0)
1506  result = result / fNormIntegral;
1507 
1508  return result;
1509  }
1510 
1511 #ifdef R__HAS_VECCORE
1512  if (fType == EFType::kTemplVec) {
1513  if (fFunctor) {
1514  if (params) result = EvalParVec(x, params);
1515  else result = EvalParVec(x, (Double_t *) fParams->GetParameters());
1516  }
1517  else {
1518  result = GetSave(x);
1519  }
1520 
1521  if (fNormalized && fNormIntegral != 0)
1522  result = result / fNormIntegral;
1523 
1524  return result;
1525  }
1526 #endif
1527 
1528  if (fType == EFType::kCompositionFcn) {
1529  if (!fComposition)
1530  Error("EvalPar", "Composition function not found");
1531 
1532  result = (*fComposition)(x, params);
1533  }
1534 
1535  return result;
1536 }
1537 
1538 ////////////////////////////////////////////////////////////////////////////////
1539 /// Execute action corresponding to one event.
1540 ///
1541 /// This member function is called when a F1 is clicked with the locator
1542 
1543 void TF1::ExecuteEvent(Int_t event, Int_t px, Int_t py)
1544 {
1545  if (!gPad) return;
1546 
1547  if (fHistogram) fHistogram->ExecuteEvent(event, px, py);
1548 
1549  if (!gPad->GetView()) {
1550  if (event == kMouseMotion) gPad->SetCursor(kHand);
1551  }
1552 }
1553 
1554 
1555 ////////////////////////////////////////////////////////////////////////////////
1556 /// Fix the value of a parameter
1557 /// The specified value will be used in a fit operation
1558 
1560 {
1561  if (ipar < 0 || ipar > GetNpar() - 1) return;
1562  SetParameter(ipar, value);
1563  if (value != 0) SetParLimits(ipar, value, value);
1564  else SetParLimits(ipar, 1, 1);
1565 }
1566 
1567 
1568 ////////////////////////////////////////////////////////////////////////////////
1569 /// Static function returning the current function being processed
1570 
1572 {
1573  ::Warning("TF1::GetCurrent", "This function is obsolete and is working only for the current painted functions");
1574  return fgCurrent;
1575 }
1576 
1577 
1578 ////////////////////////////////////////////////////////////////////////////////
1579 /// Return a pointer to the histogram used to visualise the function
1580 
1582 {
1583  if (fHistogram) return fHistogram;
1584 
1585  // histogram has not been yet created - create it
1586  // should not we make this function not const ??
1587  const_cast<TF1 *>(this)->fHistogram = const_cast<TF1 *>(this)->CreateHistogram();
1588  if (!fHistogram) Error("GetHistogram", "Error creating histogram for function %s of type %s", GetName(), IsA()->GetName());
1589  return fHistogram;
1590 }
1591 
1592 
1593 ////////////////////////////////////////////////////////////////////////////////
1594 /// Returns the maximum value of the function
1595 ///
1596 /// Method:
1597 /// First, the grid search is used to bracket the maximum
1598 /// with the step size = (xmax-xmin)/fNpx.
1599 /// This way, the step size can be controlled via the SetNpx() function.
1600 /// If the function is unimodal or if its extrema are far apart, setting
1601 /// the fNpx to a small value speeds the algorithm up many times.
1602 /// Then, Brent's method is applied on the bracketed interval
1603 /// epsilon (default = 1.E-10) controls the relative accuracy (if |x| > 1 )
1604 /// and absolute (if |x| < 1) and maxiter (default = 100) controls the maximum number
1605 /// of iteration of the Brent algorithm
1606 /// If the flag logx is set the grid search is done in log step size
1607 /// This is done automatically if the log scale is set in the current Pad
1608 ///
1609 /// NOTE: see also TF1::GetMaximumX and TF1::GetX
1610 
1612 {
1613  if (xmin >= xmax) {
1614  xmin = fXmin;
1615  xmax = fXmax;
1616  }
1617 
1618  if (!logx && gPad != 0) logx = gPad->GetLogx();
1619 
1621  GInverseFunc g(this);
1623  bm.SetFunction(wf1, xmin, xmax);
1624  bm.SetNpx(fNpx);
1625  bm.SetLogScan(logx);
1626  bm.Minimize(maxiter, epsilon, epsilon);
1627  Double_t x;
1628  x = - bm.FValMinimum();
1629 
1630  return x;
1631 }
1632 
1633 
1634 ////////////////////////////////////////////////////////////////////////////////
1635 /// Returns the X value corresponding to the maximum value of the function
1636 ///
1637 /// Method:
1638 /// First, the grid search is used to bracket the maximum
1639 /// with the step size = (xmax-xmin)/fNpx.
1640 /// This way, the step size can be controlled via the SetNpx() function.
1641 /// If the function is unimodal or if its extrema are far apart, setting
1642 /// the fNpx to a small value speeds the algorithm up many times.
1643 /// Then, Brent's method is applied on the bracketed interval
1644 /// epsilon (default = 1.E-10) controls the relative accuracy (if |x| > 1 )
1645 /// and absolute (if |x| < 1) and maxiter (default = 100) controls the maximum number
1646 /// of iteration of the Brent algorithm
1647 /// If the flag logx is set the grid search is done in log step size
1648 /// This is done automatically if the log scale is set in the current Pad
1649 ///
1650 /// NOTE: see also TF1::GetX
1651 
1653 {
1654  if (xmin >= xmax) {
1655  xmin = fXmin;
1656  xmax = fXmax;
1657  }
1658 
1659  if (!logx && gPad != 0) logx = gPad->GetLogx();
1660 
1662  GInverseFunc g(this);
1664  bm.SetFunction(wf1, xmin, xmax);
1665  bm.SetNpx(fNpx);
1666  bm.SetLogScan(logx);
1667  bm.Minimize(maxiter, epsilon, epsilon);
1668  Double_t x;
1669  x = bm.XMinimum();
1670 
1671  return x;
1672 }
1673 
1674 
1675 ////////////////////////////////////////////////////////////////////////////////
1676 /// Returns the minimum value of the function on the (xmin, xmax) interval
1677 ///
1678 /// Method:
1679 /// First, the grid search is used to bracket the maximum
1680 /// with the step size = (xmax-xmin)/fNpx. This way, the step size
1681 /// can be controlled via the SetNpx() function. If the function is
1682 /// unimodal or if its extrema are far apart, setting the fNpx to
1683 /// a small value speeds the algorithm up many times.
1684 /// Then, Brent's method is applied on the bracketed interval
1685 /// epsilon (default = 1.E-10) controls the relative accuracy (if |x| > 1 )
1686 /// and absolute (if |x| < 1) and maxiter (default = 100) controls the maximum number
1687 /// of iteration of the Brent algorithm
1688 /// If the flag logx is set the grid search is done in log step size
1689 /// This is done automatically if the log scale is set in the current Pad
1690 ///
1691 /// NOTE: see also TF1::GetMaximumX and TF1::GetX
1692 
1694 {
1695  if (xmin >= xmax) {
1696  xmin = fXmin;
1697  xmax = fXmax;
1698  }
1699 
1700  if (!logx && gPad != 0) logx = gPad->GetLogx();
1701 
1704  bm.SetFunction(wf1, xmin, xmax);
1705  bm.SetNpx(fNpx);
1706  bm.SetLogScan(logx);
1707  bm.Minimize(maxiter, epsilon, epsilon);
1708  Double_t x;
1709  x = bm.FValMinimum();
1710 
1711  return x;
1712 }
1713 
1714 ////////////////////////////////////////////////////////////////////////////////
1715 /// Find the minimum of a function of whatever dimension.
1716 /// While GetMinimum works only for 1D function , GetMinimumNDim works for all dimensions
1717 /// since it uses the minimizer interface
1718 /// vector x at beginning will contained the initial point, on exit will contain the result
1719 
1720 Double_t TF1::GetMinMaxNDim(Double_t *x , bool findmax, Double_t epsilon, Int_t maxiter) const
1721 {
1722  R__ASSERT(x != 0);
1723 
1724  int ndim = GetNdim();
1725  if (ndim == 0) {
1726  Error("GetMinimumNDim", "Function of dimension 0 - return Eval(x)");
1727  return (const_cast<TF1 &>(*this))(x);
1728  }
1729 
1730  // create minimizer class
1731  const char *minimName = ROOT::Math::MinimizerOptions::DefaultMinimizerType().c_str();
1732  const char *minimAlgo = ROOT::Math::MinimizerOptions::DefaultMinimizerAlgo().c_str();
1733  ROOT::Math::Minimizer *min = ROOT::Math::Factory::CreateMinimizer(minimName, minimAlgo);
1734 
1735  if (min == 0) {
1736  Error("GetMinimumNDim", "Error creating minimizer %s", minimName);
1737  return 0;
1738  }
1739 
1740  // minimizer will be set using default values
1741  if (epsilon > 0) min->SetTolerance(epsilon);
1742  if (maxiter > 0) min->SetMaxFunctionCalls(maxiter);
1743 
1744  // create wrapper class from TF1 (cannot use Functor, t.b.i.)
1745  ROOT::Math::WrappedMultiFunction<TF1 &> objFunc(const_cast<TF1 &>(*this), ndim);
1746  // create -f(x) when searching for the maximum
1747  GInverseFuncNdim invFunc(const_cast<TF1 *>(this));
1748  ROOT::Math::WrappedMultiFunction<GInverseFuncNdim &> objFuncInv(invFunc, ndim);
1749  if (!findmax)
1750  min->SetFunction(objFunc);
1751  else
1752  min->SetFunction(objFuncInv);
1753 
1754  std::vector<double> rmin(ndim);
1755  std::vector<double> rmax(ndim);
1756  GetRange(&rmin[0], &rmax[0]);
1757  for (int i = 0; i < ndim; ++i) {
1758  const char *xname = 0;
1759  double stepSize = 0.1;
1760  // use range for step size or give some value depending on x if range is not defined
1761  if (rmax[i] > rmin[i])
1762  stepSize = (rmax[i] - rmin[i]) / 100;
1763  else if (std::abs(x[i]) > 1.)
1764  stepSize = 0.1 * x[i];
1765 
1766  // set variable names
1767  if (ndim <= 3) {
1768  if (i == 0) {
1769  xname = "x";
1770  } else if (i == 1) {
1771  xname = "y";
1772  } else {
1773  xname = "z";
1774  }
1775  } else {
1776  xname = TString::Format("x_%d", i);
1777  // arbitrary step sie (should be computed from range)
1778  }
1779 
1780  if (rmin[i] < rmax[i]) {
1781  //Info("GetMinMax","setting limits on %s - [ %f , %f ]",xname,rmin[i],rmax[i]);
1782  min->SetLimitedVariable(i, xname, x[i], stepSize, rmin[i], rmax[i]);
1783  } else {
1784  min->SetVariable(i, xname, x[i], stepSize);
1785  }
1786  }
1787 
1788  bool ret = min->Minimize();
1789  if (!ret) {
1790  Error("GetMinimumNDim", "Error minimizing function %s", GetName());
1791  }
1792  if (min->X()) std::copy(min->X(), min->X() + ndim, x);
1793  double fmin = min->MinValue();
1794  delete min;
1795  // need to revert sign in case looking for maximum
1796  return (findmax) ? -fmin : fmin;
1797 
1798 }
1799 
1800 
1801 ////////////////////////////////////////////////////////////////////////////////
1802 /// Returns the X value corresponding to the minimum value of the function
1803 /// on the (xmin, xmax) interval
1804 ///
1805 /// Method:
1806 /// First, the grid search is used to bracket the maximum
1807 /// with the step size = (xmax-xmin)/fNpx. This way, the step size
1808 /// can be controlled via the SetNpx() function. If the function is
1809 /// unimodal or if its extrema are far apart, setting the fNpx to
1810 /// a small value speeds the algorithm up many times.
1811 /// Then, Brent's method is applied on the bracketed interval
1812 /// epsilon (default = 1.E-10) controls the relative accuracy (if |x| > 1 )
1813 /// and absolute (if |x| < 1) and maxiter (default = 100) controls the maximum number
1814 /// of iteration of the Brent algorithm
1815 /// If the flag logx is set the grid search is done in log step size
1816 /// This is done automatically if the log scale is set in the current Pad
1817 ///
1818 /// NOTE: see also TF1::GetX
1819 
1821 {
1822  if (xmin >= xmax) {
1823  xmin = fXmin;
1824  xmax = fXmax;
1825  }
1826 
1829  bm.SetFunction(wf1, xmin, xmax);
1830  bm.SetNpx(fNpx);
1831  bm.SetLogScan(logx);
1832  bm.Minimize(maxiter, epsilon, epsilon);
1833  Double_t x;
1834  x = bm.XMinimum();
1835 
1836  return x;
1837 }
1838 
1839 
1840 ////////////////////////////////////////////////////////////////////////////////
1841 /// Returns the X value corresponding to the function value fy for (xmin<x<xmax).
1842 /// in other words it can find the roots of the function when fy=0 and successive calls
1843 /// by changing the next call to [xmin+eps,xmax] where xmin is the previous root.
1844 ///
1845 /// Method:
1846 /// First, the grid search is used to bracket the maximum
1847 /// with the step size = (xmax-xmin)/fNpx. This way, the step size
1848 /// can be controlled via the SetNpx() function. If the function is
1849 /// unimodal or if its extrema are far apart, setting the fNpx to
1850 /// a small value speeds the algorithm up many times.
1851 /// Then, Brent's method is applied on the bracketed interval
1852 /// epsilon (default = 1.E-10) controls the relative accuracy (if |x| > 1 )
1853 /// and absolute (if |x| < 1) and maxiter (default = 100) controls the maximum number
1854 /// of iteration of the Brent algorithm
1855 /// If the flag logx is set the grid search is done in log step size
1856 /// This is done automatically if the log scale is set in the current Pad
1857 ///
1858 /// NOTE: see also TF1::GetMaximumX, TF1::GetMinimumX
1859 
1861 {
1862  if (xmin >= xmax) {
1863  xmin = fXmin;
1864  xmax = fXmax;
1865  }
1866 
1867  if (!logx && gPad != 0) logx = gPad->GetLogx();
1868 
1869  GFunc g(this, fy);
1872  brf.SetFunction(wf1, xmin, xmax);
1873  brf.SetNpx(fNpx);
1874  brf.SetLogScan(logx);
1875  bool ret = brf.Solve(maxiter, epsilon, epsilon);
1876  if (!ret) Error("GetX","[%f,%f] is not a valid interval",xmin,xmax);
1877  return (ret) ? brf.Root() : TMath::QuietNaN();
1878 }
1879 
1880 ////////////////////////////////////////////////////////////////////////////////
1881 /// Return the number of degrees of freedom in the fit
1882 /// the fNDF parameter has been previously computed during a fit.
1883 /// The number of degrees of freedom corresponds to the number of points
1884 /// used in the fit minus the number of free parameters.
1885 
1887 {
1888  Int_t npar = GetNpar();
1889  if (fNDF == 0 && (fNpfits > npar)) return fNpfits - npar;
1890  return fNDF;
1891 }
1892 
1893 
1894 ////////////////////////////////////////////////////////////////////////////////
1895 /// Return the number of free parameters
1896 
1898 {
1899  Int_t ntot = GetNpar();
1900  Int_t nfree = ntot;
1901  Double_t al, bl;
1902  for (Int_t i = 0; i < ntot; i++) {
1903  ((TF1 *)this)->GetParLimits(i, al, bl);
1904  if (al * bl != 0 && al >= bl) nfree--;
1905  }
1906  return nfree;
1907 }
1908 
1909 
1910 ////////////////////////////////////////////////////////////////////////////////
1911 /// Redefines TObject::GetObjectInfo.
1912 /// Displays the function info (x, function value)
1913 /// corresponding to cursor position px,py
1914 
1915 char *TF1::GetObjectInfo(Int_t px, Int_t /* py */) const
1916 {
1917  static char info[64];
1918  Double_t x = gPad->PadtoX(gPad->AbsPixeltoX(px));
1919  snprintf(info, 64, "(x=%g, f=%g)", x, ((TF1 *)this)->Eval(x));
1920  return info;
1921 }
1922 
1923 
1924 ////////////////////////////////////////////////////////////////////////////////
1925 /// Return value of parameter number ipar
1926 
1928 {
1929  if (ipar < 0 || ipar > GetNpar() - 1) return 0;
1930  return fParErrors[ipar];
1931 }
1932 
1933 
1934 ////////////////////////////////////////////////////////////////////////////////
1935 /// Return limits for parameter ipar.
1936 
1937 void TF1::GetParLimits(Int_t ipar, Double_t &parmin, Double_t &parmax) const
1938 {
1939  parmin = 0;
1940  parmax = 0;
1941  int n = fParMin.size();
1942  assert(n == int(fParMax.size()) && n <= fNpar);
1943  if (ipar < 0 || ipar > n - 1) return;
1944  parmin = fParMin[ipar];
1945  parmax = fParMax[ipar];
1946 }
1947 
1948 
1949 ////////////////////////////////////////////////////////////////////////////////
1950 /// Return the fit probability
1951 
1953 {
1954  if (fNDF <= 0) return 0;
1955  return TMath::Prob(fChisquare, fNDF);
1956 }
1957 
1958 
1959 ////////////////////////////////////////////////////////////////////////////////
1960 /// Compute Quantiles for density distribution of this function
1961 ///
1962 /// Quantile x_q of a probability distribution Function F is defined as
1963 /// \f[
1964 /// F(x_{q}) = \int_{xmin}^{x_{q}} f dx = q with 0 <= q <= 1.
1965 /// \f]
1966 /// For instance the median \f$ x_{\frac{1}{2}} \f$ of a distribution is defined as that value
1967 /// of the random variable for which the distribution function equals 0.5:
1968 /// \f[
1969 /// F(x_{\frac{1}{2}}) = \prod(x < x_{\frac{1}{2}}) = \frac{1}{2}
1970 /// \f]
1971 ///
1972 /// \param[in] this TF1 function
1973 /// \param[in] nprobSum maximum size of array q and size of array probSum
1974 /// \param[in] probSum array of positions where quantiles will be computed.
1975 /// It is assumed to contain at least nprobSum values.
1976 /// \param[out] return value nq (<=nprobSum) with the number of quantiles computed
1977 /// \param[out] array q filled with nq quantiles
1978 ///
1979 /// Getting quantiles from two histograms and storing results in a TGraph,
1980 /// a so-called QQ-plot
1981 ///
1982 /// TGraph *gr = new TGraph(nprob);
1983 /// f1->GetQuantiles(nprob,gr->GetX());
1984 /// f2->GetQuantiles(nprob,gr->GetY());
1985 /// gr->Draw("alp");
1986 ///
1987 /// \author Eddy Offermann
1988 
1989 
1990 Int_t TF1::GetQuantiles(Int_t nprobSum, Double_t *q, const Double_t *probSum)
1991 {
1992  // LM: change to use fNpx
1993  // should we change code to use a root finder ?
1994  // It should be more precise and more efficient
1995  const Int_t npx = TMath::Max(fNpx, 2 * nprobSum);
1996  const Double_t xMin = GetXmin();
1997  const Double_t xMax = GetXmax();
1998  const Double_t dx = (xMax - xMin) / npx;
1999 
2000  TArrayD integral(npx + 1);
2001  TArrayD alpha(npx);
2002  TArrayD beta(npx);
2003  TArrayD gamma(npx);
2004 
2005  integral[0] = 0;
2006  Int_t intNegative = 0;
2007  Int_t i;
2008  for (i = 0; i < npx; i++) {
2009  Double_t integ = Integral(Double_t(xMin + i * dx), Double_t(xMin + i * dx + dx), 0.0);
2010  if (integ < 0) {
2011  intNegative++;
2012  integ = -integ;
2013  }
2014  integral[i + 1] = integral[i] + integ;
2015  }
2016 
2017  if (intNegative > 0)
2018  Warning("GetQuantiles", "function:%s has %d negative values: abs assumed",
2019  GetName(), intNegative);
2020  if (integral[npx] == 0) {
2021  Error("GetQuantiles", "Integral of function is zero");
2022  return 0;
2023  }
2024 
2025  const Double_t total = integral[npx];
2026  for (i = 1; i <= npx; i++) integral[i] /= total;
2027  //the integral r for each bin is approximated by a parabola
2028  // x = alpha + beta*r +gamma*r**2
2029  // compute the coefficients alpha, beta, gamma for each bin
2030  for (i = 0; i < npx; i++) {
2031  const Double_t x0 = xMin + dx * i;
2032  const Double_t r2 = integral[i + 1] - integral[i];
2033  const Double_t r1 = Integral(x0, x0 + 0.5 * dx, 0.0) / total;
2034  gamma[i] = (2 * r2 - 4 * r1) / (dx * dx);
2035  beta[i] = r2 / dx - gamma[i] * dx;
2036  alpha[i] = x0;
2037  gamma[i] *= 2;
2038  }
2039 
2040  // Be careful because of finite precision in the integral; Use the fact that the integral
2041  // is monotone increasing
2042  for (i = 0; i < nprobSum; i++) {
2043  const Double_t r = probSum[i];
2044  Int_t bin = TMath::Max(TMath::BinarySearch(npx + 1, integral.GetArray(), r), (Long64_t)0);
2045  // in case the prob is 1
2046  if (bin == npx) {
2047  q[i] = xMax;
2048  continue;
2049  }
2050  // LM use a tolerance 1.E-12 (integral precision)
2051  while (bin < npx - 1 && TMath::AreEqualRel(integral[bin + 1], r, 1E-12)) {
2052  if (TMath::AreEqualRel(integral[bin + 2], r, 1E-12)) bin++;
2053  else break;
2054  }
2055 
2056  const Double_t rr = r - integral[bin];
2057  if (rr != 0.0) {
2058  Double_t xx = 0.0;
2059  const Double_t fac = -2.*gamma[bin] * rr / beta[bin] / beta[bin];
2060  if (fac != 0 && fac <= 1)
2061  xx = (-beta[bin] + TMath::Sqrt(beta[bin] * beta[bin] + 2 * gamma[bin] * rr)) / gamma[bin];
2062  else if (beta[bin] != 0.)
2063  xx = rr / beta[bin];
2064  q[i] = alpha[bin] + xx;
2065  } else {
2066  q[i] = alpha[bin];
2067  if (integral[bin + 1] == r) q[i] += dx;
2068  }
2069  }
2070 
2071  return nprobSum;
2072 }
2073 
2074 
2075 ////////////////////////////////////////////////////////////////////////////////
2076 /// Return a random number following this function shape
2077 ///
2078 /// The distribution contained in the function fname (TF1) is integrated
2079 /// over the channel contents.
2080 /// It is normalized to 1.
2081 /// For each bin the integral is approximated by a parabola.
2082 /// The parabola coefficients are stored as non persistent data members
2083 /// Getting one random number implies:
2084 /// - Generating a random number between 0 and 1 (say r1)
2085 /// - Look in which bin in the normalized integral r1 corresponds to
2086 /// - Evaluate the parabolic curve in the selected bin to find the corresponding X value.
2087 ///
2088 /// If the ratio fXmax/fXmin > fNpx the integral is tabulated in log scale in x
2089 /// The parabolic approximation is very good as soon as the number of bins is greater than 50.
2090 
2092 {
2093  // Check if integral array must be build
2094  if (fIntegral.size() == 0) {
2095  fIntegral.resize(fNpx + 1);
2096  fAlpha.resize(fNpx + 1);
2097  fBeta.resize(fNpx);
2098  fGamma.resize(fNpx);
2099  fIntegral[0] = 0;
2100  fAlpha[fNpx] = 0;
2101  Double_t integ;
2102  Int_t intNegative = 0;
2103  Int_t i;
2104  Bool_t logbin = kFALSE;
2105  Double_t dx;
2106  Double_t xmin = fXmin;
2107  Double_t xmax = fXmax;
2108  if (xmin > 0 && xmax / xmin > fNpx) {
2109  logbin = kTRUE;
2110  fAlpha[fNpx] = 1;
2111  xmin = TMath::Log10(fXmin);
2112  xmax = TMath::Log10(fXmax);
2113  }
2114  dx = (xmax - xmin) / fNpx;
2115 
2116  Double_t *xx = new Double_t[fNpx + 1];
2117  for (i = 0; i < fNpx; i++) {
2118  xx[i] = xmin + i * dx;
2119  }
2120  xx[fNpx] = xmax;
2121  for (i = 0; i < fNpx; i++) {
2122  if (logbin) {
2123  integ = Integral(TMath::Power(10, xx[i]), TMath::Power(10, xx[i + 1]), 0.0);
2124  } else {
2125  integ = Integral(xx[i], xx[i + 1], 0.0);
2126  }
2127  if (integ < 0) {
2128  intNegative++;
2129  integ = -integ;
2130  }
2131  fIntegral[i + 1] = fIntegral[i] + integ;
2132  }
2133  if (intNegative > 0) {
2134  Warning("GetRandom", "function:%s has %d negative values: abs assumed", GetName(), intNegative);
2135  }
2136  if (fIntegral[fNpx] == 0) {
2137  delete [] xx;
2138  Error("GetRandom", "Integral of function is zero");
2139  return 0;
2140  }
2142  for (i = 1; i <= fNpx; i++) { // normalize integral to 1
2143  fIntegral[i] /= total;
2144  }
2145  //the integral r for each bin is approximated by a parabola
2146  // x = alpha + beta*r +gamma*r**2
2147  // compute the coefficients alpha, beta, gamma for each bin
2148  Double_t x0, r1, r2, r3;
2149  for (i = 0; i < fNpx; i++) {
2150  x0 = xx[i];
2151  r2 = fIntegral[i + 1] - fIntegral[i];
2152  if (logbin) r1 = Integral(TMath::Power(10, x0), TMath::Power(10, x0 + 0.5 * dx), 0.0) / total;
2153  else r1 = Integral(x0, x0 + 0.5 * dx, 0.0) / total;
2154  r3 = 2 * r2 - 4 * r1;
2155  if (TMath::Abs(r3) > 1e-8) fGamma[i] = r3 / (dx * dx);
2156  else fGamma[i] = 0;
2157  fBeta[i] = r2 / dx - fGamma[i] * dx;
2158  fAlpha[i] = x0;
2159  fGamma[i] *= 2;
2160  }
2161  delete [] xx;
2162  }
2163 
2164  // return random number
2165  Double_t r = gRandom->Rndm();
2166  Int_t bin = TMath::BinarySearch(fNpx, fIntegral.data(), r);
2167  Double_t rr = r - fIntegral[bin];
2168 
2169  Double_t yy;
2170  if (fGamma[bin] != 0)
2171  yy = (-fBeta[bin] + TMath::Sqrt(fBeta[bin] * fBeta[bin] + 2 * fGamma[bin] * rr)) / fGamma[bin];
2172  else
2173  yy = rr / fBeta[bin];
2174  Double_t x = fAlpha[bin] + yy;
2175  if (fAlpha[fNpx] > 0) return TMath::Power(10, x);
2176  return x;
2177 }
2178 
2179 
2180 ////////////////////////////////////////////////////////////////////////////////
2181 /// Return a random number following this function shape in [xmin,xmax]
2182 ///
2183 /// The distribution contained in the function fname (TF1) is integrated
2184 /// over the channel contents.
2185 /// It is normalized to 1.
2186 /// For each bin the integral is approximated by a parabola.
2187 /// The parabola coefficients are stored as non persistent data members
2188 /// Getting one random number implies:
2189 /// - Generating a random number between 0 and 1 (say r1)
2190 /// - Look in which bin in the normalized integral r1 corresponds to
2191 /// - Evaluate the parabolic curve in the selected bin to find
2192 /// the corresponding X value.
2193 ///
2194 /// The parabolic approximation is very good as soon as the number
2195 /// of bins is greater than 50.
2196 ///
2197 /// IMPORTANT NOTE
2198 ///
2199 /// The integral of the function is computed at fNpx points. If the function
2200 /// has sharp peaks, you should increase the number of points (SetNpx)
2201 /// such that the peak is correctly tabulated at several points.
2202 
2204 {
2205  // Check if integral array must be build
2206  if (fIntegral.size() == 0) {
2207  fIntegral.resize(fNpx + 1);
2208  fAlpha.resize(fNpx+1);
2209  fBeta.resize(fNpx);
2210  fGamma.resize(fNpx);
2211 
2212  Double_t dx = (fXmax - fXmin) / fNpx;
2213  Double_t integ;
2214  Int_t intNegative = 0;
2215  Int_t i;
2216  for (i = 0; i < fNpx; i++) {
2217  integ = Integral(Double_t(fXmin + i * dx), Double_t(fXmin + i * dx + dx), 0.0);
2218  if (integ < 0) {
2219  intNegative++;
2220  integ = -integ;
2221  }
2222  fIntegral[i + 1] = fIntegral[i] + integ;
2223  }
2224  if (intNegative > 0) {
2225  Warning("GetRandom", "function:%s has %d negative values: abs assumed", GetName(), intNegative);
2226  }
2227  if (fIntegral[fNpx] == 0) {
2228  Error("GetRandom", "Integral of function is zero");
2229  return 0;
2230  }
2232  for (i = 1; i <= fNpx; i++) { // normalize integral to 1
2233  fIntegral[i] /= total;
2234  }
2235  //the integral r for each bin is approximated by a parabola
2236  // x = alpha + beta*r +gamma*r**2
2237  // compute the coefficients alpha, beta, gamma for each bin
2238  Double_t x0, r1, r2, r3;
2239  for (i = 0; i < fNpx; i++) {
2240  x0 = fXmin + i * dx;
2241  r2 = fIntegral[i + 1] - fIntegral[i];
2242  r1 = Integral(x0, x0 + 0.5 * dx, 0.0) / total;
2243  r3 = 2 * r2 - 4 * r1;
2244  if (TMath::Abs(r3) > 1e-8) fGamma[i] = r3 / (dx * dx);
2245  else fGamma[i] = 0;
2246  fBeta[i] = r2 / dx - fGamma[i] * dx;
2247  fAlpha[i] = x0;
2248  fGamma[i] *= 2;
2249  }
2250  }
2251 
2252  // return random number
2253  Double_t dx = (fXmax - fXmin) / fNpx;
2254  Int_t nbinmin = (Int_t)((xmin - fXmin) / dx);
2255  Int_t nbinmax = (Int_t)((xmax - fXmin) / dx) + 2;
2256  if (nbinmax > fNpx) nbinmax = fNpx;
2257 
2258  Double_t pmin = fIntegral[nbinmin];
2259  Double_t pmax = fIntegral[nbinmax];
2260 
2261  Double_t r, x, xx, rr;
2262  do {
2263  r = gRandom->Uniform(pmin, pmax);
2264 
2265  Int_t bin = TMath::BinarySearch(fNpx, fIntegral.data(), r);
2266  rr = r - fIntegral[bin];
2267 
2268  if (fGamma[bin] != 0)
2269  xx = (-fBeta[bin] + TMath::Sqrt(fBeta[bin] * fBeta[bin] + 2 * fGamma[bin] * rr)) / fGamma[bin];
2270  else
2271  xx = rr / fBeta[bin];
2272  x = fAlpha[bin] + xx;
2273  } while (x < xmin || x > xmax);
2274  return x;
2275 }
2276 
2277 ////////////////////////////////////////////////////////////////////////////////
2278 /// Return range of a generic N-D function.
2279 
2280 void TF1::GetRange(Double_t *rmin, Double_t *rmax) const
2281 {
2282  int ndim = GetNdim();
2283 
2284  double xmin = 0, ymin = 0, zmin = 0, xmax = 0, ymax = 0, zmax = 0;
2285  GetRange(xmin, ymin, zmin, xmax, ymax, zmax);
2286  for (int i = 0; i < ndim; ++i) {
2287  if (i == 0) {
2288  rmin[0] = xmin;
2289  rmax[0] = xmax;
2290  } else if (i == 1) {
2291  rmin[1] = ymin;
2292  rmax[1] = ymax;
2293  } else if (i == 2) {
2294  rmin[2] = zmin;
2295  rmax[2] = zmax;
2296  } else {
2297  rmin[i] = 0;
2298  rmax[i] = 0;
2299  }
2300  }
2301 }
2302 
2303 
2304 ////////////////////////////////////////////////////////////////////////////////
2305 /// Return range of a 1-D function.
2306 
2308 {
2309  xmin = fXmin;
2310  xmax = fXmax;
2311 }
2312 
2313 
2314 ////////////////////////////////////////////////////////////////////////////////
2315 /// Return range of a 2-D function.
2316 
2318 {
2319  xmin = fXmin;
2320  xmax = fXmax;
2321  ymin = 0;
2322  ymax = 0;
2323 }
2324 
2325 
2326 ////////////////////////////////////////////////////////////////////////////////
2327 /// Return range of function.
2328 
2330 {
2331  xmin = fXmin;
2332  xmax = fXmax;
2333  ymin = 0;
2334  ymax = 0;
2335  zmin = 0;
2336  zmax = 0;
2337 }
2338 
2339 
2340 ////////////////////////////////////////////////////////////////////////////////
2341 /// Get value corresponding to X in array of fSave values
2342 
2344 {
2345  if (fSave.size() == 0) return 0;
2346  //if (fSave == 0) return 0;
2347  int fNsave = fSave.size();
2348  Double_t x = Double_t(xx[0]);
2349  Double_t y, dx, xmin, xmax, xlow, xup, ylow, yup;
2350  if (fParent && fParent->InheritsFrom(TH1::Class())) {
2351  //if parent is a histogram the function had been saved at the center of the bins
2352  //we make a linear interpolation between the saved values
2353  xmin = fSave[fNsave - 3];
2354  xmax = fSave[fNsave - 2];
2355  if (fSave[fNsave - 1] == xmax) {
2356  TH1 *h = (TH1 *)fParent;
2357  TAxis *xaxis = h->GetXaxis();
2358  Int_t bin1 = xaxis->FindBin(xmin);
2359  Int_t binup = xaxis->FindBin(xmax);
2360  Int_t bin = xaxis->FindBin(x);
2361  if (bin < binup) {
2362  xlow = xaxis->GetBinCenter(bin);
2363  xup = xaxis->GetBinCenter(bin + 1);
2364  ylow = fSave[bin - bin1];
2365  yup = fSave[bin - bin1 + 1];
2366  } else {
2367  xlow = xaxis->GetBinCenter(bin - 1);
2368  xup = xaxis->GetBinCenter(bin);
2369  ylow = fSave[bin - bin1 - 1];
2370  yup = fSave[bin - bin1];
2371  }
2372  dx = xup - xlow;
2373  y = ((xup * ylow - xlow * yup) + x * (yup - ylow)) / dx;
2374  return y;
2375  }
2376  }
2377  Int_t np = fNsave - 3;
2378  xmin = Double_t(fSave[np + 1]);
2379  xmax = Double_t(fSave[np + 2]);
2380  dx = (xmax - xmin) / np;
2381  if (x < xmin || x > xmax) return 0;
2382  // return a Nan in case of x=nan, otherwise will crash later
2383  if (TMath::IsNaN(x)) return x;
2384  if (dx <= 0) return 0;
2385 
2386  Int_t bin = Int_t((x - xmin) / dx);
2387  xlow = xmin + bin * dx;
2388  xup = xlow + dx;
2389  ylow = fSave[bin];
2390  yup = fSave[bin + 1];
2391  y = ((xup * ylow - xlow * yup) + x * (yup - ylow)) / dx;
2392  return y;
2393 }
2394 
2395 
2396 ////////////////////////////////////////////////////////////////////////////////
2397 /// Get x axis of the function.
2398 
2400 {
2401  TH1 *h = GetHistogram();
2402  if (!h) return 0;
2403  return h->GetXaxis();
2404 }
2405 
2406 
2407 ////////////////////////////////////////////////////////////////////////////////
2408 /// Get y axis of the function.
2409 
2411 {
2412  TH1 *h = GetHistogram();
2413  if (!h) return 0;
2414  return h->GetYaxis();
2415 }
2416 
2417 
2418 ////////////////////////////////////////////////////////////////////////////////
2419 /// Get z axis of the function. (In case this object is a TF2 or TF3)
2420 
2422 {
2423  TH1 *h = GetHistogram();
2424  if (!h) return 0;
2425  return h->GetZaxis();
2426 }
2427 
2428 
2429 
2430 ////////////////////////////////////////////////////////////////////////////////
2431 /// Compute the gradient (derivative) wrt a parameter ipar
2432 ///
2433 /// \param ipar index of parameter for which the derivative is computed
2434 /// \param x point, where the derivative is computed
2435 /// \param eps - if the errors of parameters have been computed, the step used in
2436 /// numerical differentiation is eps*parameter_error.
2437 ///
2438 /// if the errors have not been computed, step=eps is used
2439 /// default value of eps = 0.01
2440 /// Method is the same as in Derivative() function
2441 ///
2442 /// If a parameter is fixed, the gradient on this parameter = 0
2443 
2445 {
2446  return GradientParTempl<Double_t>(ipar, x, eps);
2447 }
2448 
2449 ////////////////////////////////////////////////////////////////////////////////
2450 /// Compute the gradient wrt parameters
2451 ///
2452 /// \param x point, were the gradient is computed
2453 /// \param grad used to return the computed gradient, assumed to be of at least fNpar size
2454 /// \param eps if the errors of parameters have been computed, the step used in
2455 /// numerical differentiation is eps*parameter_error.
2456 ///
2457 /// if the errors have not been computed, step=eps is used
2458 /// default value of eps = 0.01
2459 /// Method is the same as in Derivative() function
2460 ///
2461 /// If a parameter is fixed, the gradient on this parameter = 0
2462 
2463 void TF1::GradientPar(const Double_t *x, Double_t *grad, Double_t eps)
2464 {
2465  GradientParTempl<Double_t>(x, grad, eps);
2466 }
2467 
2468 ////////////////////////////////////////////////////////////////////////////////
2469 /// Initialize parameters addresses.
2470 
2471 void TF1::InitArgs(const Double_t *x, const Double_t *params)
2472 {
2473  if (fMethodCall) {
2474  Long_t args[2];
2475  args[0] = (Long_t)x;
2476  if (params) args[1] = (Long_t)params;
2477  else args[1] = (Long_t)GetParameters();
2478  fMethodCall->SetParamPtrs(args);
2479  }
2480 }
2481 
2482 
2483 ////////////////////////////////////////////////////////////////////////////////
2484 /// Create the basic function objects
2485 
2487 {
2488  TF1 *f1;
2490  if (!gROOT->GetListOfFunctions()->FindObject("gaus")) {
2491  f1 = new TF1("gaus", "gaus", -1, 1);
2492  f1->SetParameters(1, 0, 1);
2493  f1 = new TF1("gausn", "gausn", -1, 1);
2494  f1->SetParameters(1, 0, 1);
2495  f1 = new TF1("landau", "landau", -1, 1);
2496  f1->SetParameters(1, 0, 1);
2497  f1 = new TF1("landaun", "landaun", -1, 1);
2498  f1->SetParameters(1, 0, 1);
2499  f1 = new TF1("expo", "expo", -1, 1);
2500  f1->SetParameters(1, 1);
2501  for (Int_t i = 0; i < 10; i++) {
2502  f1 = new TF1(Form("pol%d", i), Form("pol%d", i), -1, 1);
2503  f1->SetParameters(1, 1, 1, 1, 1, 1, 1, 1, 1, 1);
2504  // create also chebyshev polynomial
2505  // (note polynomial object will not be deleted)
2506  // note that these functions cannot be stored
2508  Double_t min = -1;
2509  Double_t max = 1;
2510  f1 = new TF1(TString::Format("chebyshev%d", i), pol, min, max, i + 1, 1);
2511  f1->SetParameters(1, 1, 1, 1, 1, 1, 1, 1, 1, 1);
2512  }
2513 
2514  }
2515 }
2516 ////////////////////////////////////////////////////////////////////////////////
2517 /// IntegralOneDim or analytical integral
2518 
2520 {
2521  Double_t error = 0;
2522  if (GetNumber() > 0) {
2523  Double_t result = 0.;
2524  if (gDebug) {
2525  Info("computing analytical integral for function %s with number %d", GetName(), GetNumber());
2526  }
2527  result = AnalyticalIntegral(this, a, b);
2528  // if it is a formula that havent been implemented in analytical integral a NaN is return
2529  if (!TMath::IsNaN(result)) return result;
2530  if (gDebug)
2531  Warning("analytical integral not available for %s - with number %d compute numerical integral", GetName(), GetNumber());
2532  }
2533  return IntegralOneDim(a, b, epsrel, epsrel, error);
2534 }
2535 
2536 ////////////////////////////////////////////////////////////////////////////////
2537 /// Return Integral of function between a and b using the given parameter values and
2538 /// relative and absolute tolerance.
2539 ///
2540 /// The default integrator defined in ROOT::Math::IntegratorOneDimOptions::DefaultIntegrator() is used
2541 /// If ROOT contains the MathMore library the default integrator is set to be
2542 /// the adaptive ROOT::Math::GSLIntegrator (based on QUADPACK) or otherwise the
2543 /// ROOT::Math::GaussIntegrator is used
2544 /// See the reference documentation of these classes for more information about the
2545 /// integration algorithms
2546 /// To change integration algorithm just do :
2547 /// ROOT::Math::IntegratorOneDimOptions::SetDefaultIntegrator(IntegratorName);
2548 /// Valid integrator names are:
2549 /// - Gauss : for ROOT::Math::GaussIntegrator
2550 /// - GaussLegendre : for ROOT::Math::GaussLegendreIntegrator
2551 /// - Adaptive : for ROOT::Math::GSLIntegrator adaptive method (QAG)
2552 /// - AdaptiveSingular : for ROOT::Math::GSLIntegrator adaptive singular method (QAGS)
2553 /// - NonAdaptive : for ROOT::Math::GSLIntegrator non adaptive (QNG)
2554 ///
2555 /// In order to use the GSL integrators one needs to have the MathMore library installed
2556 ///
2557 /// Note 1:
2558 ///
2559 /// Values of the function f(x) at the interval end-points A and B are not
2560 /// required. The subprogram may therefore be used when these values are
2561 /// undefined.
2562 ///
2563 /// Note 2:
2564 ///
2565 /// Instead of TF1::Integral, you may want to use the combination of
2566 /// TF1::CalcGaussLegendreSamplingPoints and TF1::IntegralFast.
2567 /// See an example with the following script:
2568 ///
2569 /// ~~~ {.cpp}
2570 /// void gint() {
2571 /// TF1 *g = new TF1("g","gaus",-5,5);
2572 /// g->SetParameters(1,0,1);
2573 /// //default gaus integration method uses 6 points
2574 /// //not suitable to integrate on a large domain
2575 /// double r1 = g->Integral(0,5);
2576 /// double r2 = g->Integral(0,1000);
2577 ///
2578 /// //try with user directives computing more points
2579 /// Int_t np = 1000;
2580 /// double *x=new double[np];
2581 /// double *w=new double[np];
2582 /// g->CalcGaussLegendreSamplingPoints(np,x,w,1e-15);
2583 /// double r3 = g->IntegralFast(np,x,w,0,5);
2584 /// double r4 = g->IntegralFast(np,x,w,0,1000);
2585 /// double r5 = g->IntegralFast(np,x,w,0,10000);
2586 /// double r6 = g->IntegralFast(np,x,w,0,100000);
2587 /// printf("g->Integral(0,5) = %g\n",r1);
2588 /// printf("g->Integral(0,1000) = %g\n",r2);
2589 /// printf("g->IntegralFast(n,x,w,0,5) = %g\n",r3);
2590 /// printf("g->IntegralFast(n,x,w,0,1000) = %g\n",r4);
2591 /// printf("g->IntegralFast(n,x,w,0,10000) = %g\n",r5);
2592 /// printf("g->IntegralFast(n,x,w,0,100000)= %g\n",r6);
2593 /// delete [] x;
2594 /// delete [] w;
2595 /// }
2596 /// ~~~
2597 ///
2598 /// This example produces the following results:
2599 ///
2600 /// ~~~ {.cpp}
2601 /// g->Integral(0,5) = 1.25331
2602 /// g->Integral(0,1000) = 1.25319
2603 /// g->IntegralFast(n,x,w,0,5) = 1.25331
2604 /// g->IntegralFast(n,x,w,0,1000) = 1.25331
2605 /// g->IntegralFast(n,x,w,0,10000) = 1.25331
2606 /// g->IntegralFast(n,x,w,0,100000)= 1.253
2607 /// ~~~
2608 
2610 {
2611  //Double_t *parameters = GetParameters();
2612  TF1_EvalWrapper wf1(this, 0, fgAbsValue);
2613  Double_t result = 0;
2614  Int_t status = 0;
2618  ROOT::Math::GaussIntegrator iod(epsabs, epsrel);
2619  iod.SetFunction(wf1);
2620  if (a != - TMath::Infinity() && b != TMath::Infinity())
2621  result = iod.Integral(a, b);
2622  else if (a == - TMath::Infinity() && b != TMath::Infinity())
2623  result = iod.IntegralLow(b);
2624  else if (a != - TMath::Infinity() && b == TMath::Infinity())
2625  result = iod.IntegralUp(a);
2626  else if (a == - TMath::Infinity() && b == TMath::Infinity())
2627  result = iod.Integral();
2628  error = iod.Error();
2629  status = iod.Status();
2630  } else {
2632  if (a != - TMath::Infinity() && b != TMath::Infinity())
2633  result = iod.Integral(a, b);
2634  else if (a == - TMath::Infinity() && b != TMath::Infinity())
2635  result = iod.IntegralLow(b);
2636  else if (a != - TMath::Infinity() && b == TMath::Infinity())
2637  result = iod.IntegralUp(a);
2638  else if (a == - TMath::Infinity() && b == TMath::Infinity())
2639  result = iod.Integral();
2640  error = iod.Error();
2641  status = iod.Status();
2642  }
2643  if (status != 0) {
2645  Warning("IntegralOneDim", "Error found in integrating function %s in [%f,%f] using %s. Result = %f +/- %f - status = %d", GetName(), a, b, igName.c_str(), result, error, status);
2646  TString msg("\t\tFunction Parameters = {");
2647  for (int ipar = 0; ipar < GetNpar(); ++ipar) {
2648  msg += TString::Format(" %s = %f ", GetParName(ipar), GetParameter(ipar));
2649  if (ipar < GetNpar() - 1) msg += TString(",");
2650  else msg += TString("}");
2651  }
2652  Info("IntegralOneDim", "%s", msg.Data());
2653  }
2654  return result;
2655 }
2656 
2657 
2658 //______________________________________________________________________________
2659 // Double_t TF1::Integral(Double_t, Double_t, Double_t, Double_t, Double_t, Double_t)
2660 // {
2661 // // Return Integral of a 2d function in range [ax,bx],[ay,by]
2662 
2663 // Error("Integral","Must be called with a TF2 only");
2664 // return 0;
2665 // }
2666 
2667 
2668 // //______________________________________________________________________________
2669 // Double_t TF1::Integral(Double_t, Double_t, Double_t, Double_t, Double_t, Double_t, Double_t, Double_t)
2670 // {
2671 // // Return Integral of a 3d function in range [ax,bx],[ay,by],[az,bz]
2672 
2673 // Error("Integral","Must be called with a TF3 only");
2674 // return 0;
2675 // }
2676 
2677 ////////////////////////////////////////////////////////////////////////////////
2678 /// Return Error on Integral of a parametric function between a and b
2679 /// due to the parameter uncertainties.
2680 /// A pointer to a vector of parameter values and to the elements of the covariance matrix (covmat)
2681 /// can be optionally passed. By default (i.e. when a zero pointer is passed) the current stored
2682 /// parameter values are used to estimate the integral error together with the covariance matrix
2683 /// from the last fit (retrieved from the global fitter instance)
2684 ///
2685 /// IMPORTANT NOTE1:
2686 ///
2687 /// When no covariance matrix is passed and in the meantime a fit is done
2688 /// using another function, the routine will signal an error and it will return zero only
2689 /// when the number of fit parameter is different than the values stored in TF1 (TF1::GetNpar() ).
2690 /// In the case that npar is the same, an incorrect result is returned.
2691 ///
2692 /// IMPORTANT NOTE2:
2693 ///
2694 /// The user must pass a pointer to the elements of the full covariance matrix
2695 /// dimensioned with the right size (npar*npar), where npar is the total number of parameters (TF1::GetNpar()),
2696 /// including also the fixed parameters. When there are fixed parameters, the pointer returned from
2697 /// TVirtualFitter::GetCovarianceMatrix() cannot be used.
2698 /// One should use the TFitResult class, as shown in the example below.
2699 ///
2700 /// To get the matrix and values from an old fit do for example:
2701 /// TFitResultPtr r = histo->Fit(func, "S");
2702 /// ..... after performing other fits on the same function do
2703 ///
2704 /// func->IntegralError(x1,x2,r->GetParams(), r->GetCovarianceMatrix()->GetMatrixArray() );
2705 
2707 {
2708  Double_t x1[1];
2709  Double_t x2[1];
2710  x1[0] = a, x2[0] = b;
2711  return ROOT::TF1Helper::IntegralError(this, 1, x1, x2, params, covmat, epsilon);
2712 }
2713 
2714 ////////////////////////////////////////////////////////////////////////////////
2715 /// Return Error on Integral of a parametric function with dimension larger tan one
2716 /// between a[] and b[] due to the parameters uncertainties.
2717 /// For a TF1 with dimension larger than 1 (for example a TF2 or TF3)
2718 /// TF1::IntegralMultiple is used for the integral calculation
2719 ///
2720 /// A pointer to a vector of parameter values and to the elements of the covariance matrix (covmat) can be optionally passed.
2721 /// By default (i.e. when a zero pointer is passed) the current stored parameter values are used to estimate the integral error
2722 /// together with the covariance matrix from the last fit (retrieved from the global fitter instance).
2723 ///
2724 /// IMPORTANT NOTE1:
2725 ///
2726 /// When no covariance matrix is passed and in the meantime a fit is done
2727 /// using another function, the routine will signal an error and it will return zero only
2728 /// when the number of fit parameter is different than the values stored in TF1 (TF1::GetNpar() ).
2729 /// In the case that npar is the same, an incorrect result is returned.
2730 ///
2731 /// IMPORTANT NOTE2:
2732 ///
2733 /// The user must pass a pointer to the elements of the full covariance matrix
2734 /// dimensioned with the right size (npar*npar), where npar is the total number of parameters (TF1::GetNpar()),
2735 /// including also the fixed parameters. When there are fixed parameters, the pointer returned from
2736 /// TVirtualFitter::GetCovarianceMatrix() cannot be used.
2737 /// One should use the TFitResult class, as shown in the example below.
2738 ///
2739 /// To get the matrix and values from an old fit do for example:
2740 /// TFitResultPtr r = histo->Fit(func, "S");
2741 /// ..... after performing other fits on the same function do
2742 ///
2743 /// func->IntegralError(x1,x2,r->GetParams(), r->GetCovarianceMatrix()->GetMatrixArray() );
2744 
2745 Double_t TF1::IntegralError(Int_t n, const Double_t *a, const Double_t *b, const Double_t *params, const Double_t *covmat, Double_t epsilon)
2746 {
2747  return ROOT::TF1Helper::IntegralError(this, n, a, b, params, covmat, epsilon);
2748 }
2749 
2750 #ifdef INTHEFUTURE
2751 ////////////////////////////////////////////////////////////////////////////////
2752 /// Gauss-Legendre integral, see CalcGaussLegendreSamplingPoints
2753 
2755 {
2756  if (!g) return 0;
2757  return IntegralFast(g->GetN(), g->GetX(), g->GetY(), a, b, params);
2758 }
2759 #endif
2760 
2761 
2762 ////////////////////////////////////////////////////////////////////////////////
2763 /// Gauss-Legendre integral, see CalcGaussLegendreSamplingPoints
2764 
2766 {
2767  // Now x and w are not used!
2768 
2769  ROOT::Math::WrappedTF1 wf1(*this);
2770  if (params)
2771  wf1.SetParameters(params);
2773  gli.SetFunction(wf1);
2774  return gli.Integral(a, b);
2775 
2776 }
2777 
2778 
2779 ////////////////////////////////////////////////////////////////////////////////
2780 /// See more general prototype below.
2781 /// This interface kept for back compatibility
2782 /// It is recommended to use the other interface where one can specify also epsabs and the maximum number of
2783 /// points
2784 
2786 {
2787  Int_t nfnevl, ifail;
2789  Double_t result = IntegralMultiple(n, a, b, maxpts, epsrel, epsrel, relerr, nfnevl, ifail);
2790  if (ifail > 0) {
2791  Warning("IntegralMultiple", "failed code=%d, ", ifail);
2792  }
2793  return result;
2794 }
2795 
2796 
2797 ////////////////////////////////////////////////////////////////////////////////
2798 /// This function computes, to an attempted specified accuracy, the value of
2799 /// the integral
2800 ///
2801 /// \param[in] n Number of dimensions [2,15]
2802 /// \param[in] a,b One-dimensional arrays of length >= N . On entry A[i], and B[i],
2803 /// contain the lower and upper limits of integration, respectively.
2804 /// \param[in] maxpts Maximum number of function evaluations to be allowed.
2805 /// maxpts >= 2^n +2*n*(n+1) +1
2806 /// if maxpts<minpts, maxpts is set to 10*minpts
2807 /// \param[in] epsrel Specified relative accuracy.
2808 /// \param[in] epsabs Specified absolute accuracy.
2809 /// The integration algorithm will attempt to reach either the relative or the absolute accuracy.
2810 /// In case the maximum function called is reached the algorithm will stop earlier without having reached
2811 /// the desired accuracy
2812 ///
2813 /// \param[out] relerr Contains, on exit, an estimation of the relative accuracy of the result.
2814 /// \param[out] nfnevl number of function evaluations performed.
2815 /// \param[out] ifail
2816 /// \parblock
2817 /// 0 Normal exit. At least minpts and at most maxpts calls to the function were performed.
2818 ///
2819 /// 1 maxpts is too small for the specified accuracy eps. The result and relerr contain the values obtainable for the
2820 /// specified value of maxpts.
2821 ///
2822 /// 3 n<2 or n>15
2823 /// \endparblock
2824 ///
2825 /// Method:
2826 ///
2827 /// The default method used is the Genz-Mallik adaptive multidimensional algorithm
2828 /// using the class ROOT::Math::AdaptiveIntegratorMultiDim (see the reference documentation of the class)
2829 ///
2830 /// Other methods can be used by setting ROOT::Math::IntegratorMultiDimOptions::SetDefaultIntegrator()
2831 /// to different integrators.
2832 /// Other possible integrators are MC integrators based on the ROOT::Math::GSLMCIntegrator class
2833 /// Possible methods are : Vegas, Miser or Plain
2834 /// IN case of MC integration the accuracy is determined by the number of function calls, one should be
2835 /// careful not to use a too large value of maxpts
2836 ///
2837 
2838 Double_t TF1::IntegralMultiple(Int_t n, const Double_t *a, const Double_t *b, Int_t maxpts, Double_t epsrel, Double_t epsabs, Double_t &relerr, Int_t &nfnevl, Int_t &ifail)
2839 {
2841 
2842  double result = 0;
2846  ROOT::Math::AdaptiveIntegratorMultiDim aimd(wf1, epsabs, epsrel, maxpts);
2847  //aimd.SetMinPts(minpts); // use default minpts ( n^2 + 2 * n * (n+1) +1 )
2848  result = aimd.Integral(a, b);
2849  relerr = aimd.RelError();
2850  nfnevl = aimd.NEval();
2851  ifail = aimd.Status();
2852  } else {
2853  // use default abs tolerance = relative tolerance
2855  result = imd.Integral(a, b);
2856  relerr = (result != 0) ? imd.Error() / std::abs(result) : imd.Error();
2857  nfnevl = 0;
2858  ifail = imd.Status();
2859  }
2860 
2861 
2862  return result;
2863 }
2864 
2865 
2866 ////////////////////////////////////////////////////////////////////////////////
2867 /// Return kTRUE if the function is valid
2868 
2870 {
2871  if (fFormula) return fFormula->IsValid();
2872  if (fMethodCall) return fMethodCall->IsValid();
2873  // function built on compiled functors are always valid by definition
2874  // (checked at compiled time)
2875  // invalid is a TF1 where the functor is null pointer and has not been saved
2876  if (!fFunctor && fSave.empty()) return kFALSE;
2877  return kTRUE;
2878 }
2879 
2880 
2881 //______________________________________________________________________________
2882 
2883 
2884 void TF1::Print(Option_t *option) const
2885 {
2886  if (fType == EFType::kFormula) {
2887  printf("Formula based function: %s \n", GetName());
2888  assert(fFormula);
2889  fFormula->Print(option);
2890  } else if (fType > 0) {
2891  if (fType == EFType::kInterpreted)
2892  printf("Interpreted based function: %s(double *x, double *p). Ndim = %d, Npar = %d \n", GetName(), GetNdim(),
2893  GetNpar());
2894  else if (fType == EFType::kCompositionFcn) {
2895  printf("Composition based function: %s. Ndim = %d, Npar = %d \n", GetName(), GetNdim(), GetNpar());
2896  if (!fComposition)
2897  printf("fComposition not found!\n"); // this would be bad
2898  } else {
2899  if (fFunctor)
2900  printf("Compiled based function: %s based on a functor object. Ndim = %d, Npar = %d\n", GetName(),
2901  GetNdim(), GetNpar());
2902  else {
2903  printf("Function based on a list of points from a compiled based function: %s. Ndim = %d, Npar = %d, Npx "
2904  "= %zu\n",
2905  GetName(), GetNdim(), GetNpar(), fSave.size());
2906  if (fSave.empty())
2907  Warning("Print", "Function %s is based on a list of points but list is empty", GetName());
2908  }
2909  }
2910  TString opt(option);
2911  opt.ToUpper();
2912  if (opt.Contains("V")) {
2913  // print list of parameters
2914  if (fNpar > 0) {
2915  printf("List of Parameters: \n");
2916  for (int i = 0; i < fNpar; ++i)
2917  printf(" %20s = %10f \n", GetParName(i), GetParameter(i));
2918  }
2919  if (!fSave.empty()) {
2920  // print list of saved points
2921  printf("List of Saved points (N=%d): \n", int(fSave.size()));
2922  for (auto &x : fSave)
2923  printf("( %10f ) ", x);
2924  printf("\n");
2925  }
2926  }
2927  }
2928  if (fHistogram) {
2929  printf("Contained histogram\n");
2930  fHistogram->Print(option);
2931  }
2932 }
2933 
2934 ////////////////////////////////////////////////////////////////////////////////
2935 /// Paint this function with its current attributes.
2936 /// The function is going to be converted in an histogram and the corresponding
2937 /// histogram is painted.
2938 /// The painted histogram can be retrieved calling afterwards the method TF1::GetHistogram()
2939 
2940 void TF1::Paint(Option_t *choptin)
2941 {
2942  fgCurrent = this;
2943 
2944  char option[32];
2945  strlcpy(option,choptin,32);
2946 
2947  TString opt = option;
2948  opt.ToLower();
2949 
2950  Bool_t optSAME = kFALSE;
2951  if (opt.Contains("same")) {
2952  opt.ReplaceAll("same","");
2953  optSAME = kTRUE;
2954  }
2955  opt.ReplaceAll(' ', "");
2956 
2957  Double_t xmin = fXmin, xmax = fXmax, pmin = fXmin, pmax = fXmax;
2958  if (gPad) {
2959  pmin = gPad->PadtoX(gPad->GetUxmin());
2960  pmax = gPad->PadtoX(gPad->GetUxmax());
2961  }
2962  if (optSAME) {
2963  if (xmax < pmin) return; // Completely outside.
2964  if (xmin > pmax) return;
2965  if (xmin < pmin) xmin = pmin;
2966  if (xmax > pmax) xmax = pmax;
2967  }
2968 
2969  // create an histogram using the function content (re-use it if already existing)
2971 
2972  char *l1 = strstr(option,"PFC"); // Automatic Fill Color
2973  char *l2 = strstr(option,"PLC"); // Automatic Line Color
2974  char *l3 = strstr(option,"PMC"); // Automatic Marker Color
2975  if (l1 || l2 || l3) {
2976  Int_t i = gPad->NextPaletteColor();
2977  if (l1) {memcpy(l1," ",3); fHistogram->SetFillColor(i);}
2978  if (l2) {memcpy(l2," ",3); fHistogram->SetLineColor(i);}
2979  if (l3) {memcpy(l3," ",3); fHistogram->SetMarkerColor(i);}
2980  }
2981 
2982  // set the optimal minimum and maximum
2983  Double_t minimum = fHistogram->GetMinimumStored();
2984  Double_t maximum = fHistogram->GetMaximumStored();
2985  if (minimum <= 0 && gPad && gPad->GetLogy()) minimum = -1111; // This can happen when switching from lin to log scale.
2986  if (gPad && gPad->GetUymin() < fHistogram->GetMinimum() &&
2987  !fHistogram->TestBit(TH1::kIsZoomed)) minimum = -1111; // This can happen after unzooming a fit.
2988  if (minimum == -1111) { // This can happen after unzooming.
2990  minimum = fHistogram->GetYaxis()->GetXmin();
2991  } else {
2992  minimum = fMinimum;
2993  // Optimize the computation of the scale in Y in case the min/max of the
2994  // function oscillate around a constant value
2995  if (minimum == -1111) {
2996  Double_t hmin;
2997  if (optSAME && gPad) hmin = gPad->GetUymin();
2998  else hmin = fHistogram->GetMinimum();
2999  if (hmin > 0) {
3000  Double_t hmax;
3001  Double_t hminpos = hmin;
3002  if (optSAME && gPad) hmax = gPad->GetUymax();
3003  else hmax = fHistogram->GetMaximum();
3004  hmin -= 0.05 * (hmax - hmin);
3005  if (hmin < 0) hmin = 0;
3006  if (hmin <= 0 && gPad && gPad->GetLogy()) hmin = hminpos;
3007  minimum = hmin;
3008  }
3009  }
3010  }
3011  fHistogram->SetMinimum(minimum);
3012  }
3013  if (maximum == -1111) {
3015  maximum = fHistogram->GetYaxis()->GetXmax();
3016  } else {
3017  maximum = fMaximum;
3018  }
3019  fHistogram->SetMaximum(maximum);
3020  }
3021 
3022 
3023  // Draw the histogram.
3024  if (!gPad) return;
3025  if (opt.Length() == 0) {
3026  if (optSAME) fHistogram->Paint("lfsame");
3027  else fHistogram->Paint("lf");
3028  } else {
3029  fHistogram->Paint(option);
3030  }
3031 }
3032 
3033 ////////////////////////////////////////////////////////////////////////////////
3034 /// Create histogram with bin content equal to function value
3035 /// computed at the bin center
3036 /// This histogram will be used to paint the function
3037 /// A re-creation is forced and a new histogram is done if recreate=true
3038 
3040 {
3041  Int_t i;
3042  Double_t xv[1];
3043 
3044  TH1 *histogram = 0;
3045 
3046 
3047  // Create a temporary histogram and fill each channel with the function value
3048  // Preserve axis titles
3049  TString xtitle = "";
3050  TString ytitle = "";
3051  char *semicol = (char *)strstr(GetTitle(), ";");
3052  if (semicol) {
3053  Int_t nxt = strlen(semicol);
3054  char *ctemp = new char[nxt];
3055  strlcpy(ctemp, semicol + 1, nxt);
3056  semicol = (char *)strstr(ctemp, ";");
3057  if (semicol) {
3058  *semicol = 0;
3059  ytitle = semicol + 1;
3060  }
3061  xtitle = ctemp;
3062  delete [] ctemp;
3063  }
3064  if (fHistogram) {
3065  // delete previous histograms if were done if done in different mode
3066  xtitle = fHistogram->GetXaxis()->GetTitle();
3067  ytitle = fHistogram->GetYaxis()->GetTitle();
3068  if (!gPad->GetLogx() && fHistogram->TestBit(TH1::kLogX)) {
3069  delete fHistogram;
3070  fHistogram = 0;
3071  recreate = kTRUE;
3072  }
3073  if (gPad->GetLogx() && !fHistogram->TestBit(TH1::kLogX)) {
3074  delete fHistogram;
3075  fHistogram = 0;
3076  recreate = kTRUE;
3077  }
3078  }
3079 
3080  if (fHistogram && !recreate) {
3081  histogram = fHistogram;
3083  } else {
3084  // If logx, we must bin in logx and not in x
3085  // otherwise in case of several decades, one gets wrong results.
3086  if (xmin > 0 && gPad && gPad->GetLogx()) {
3087  Double_t *xbins = new Double_t[fNpx + 1];
3088  Double_t xlogmin = TMath::Log10(xmin);
3089  Double_t xlogmax = TMath::Log10(xmax);
3090  Double_t dlogx = (xlogmax - xlogmin) / ((Double_t)fNpx);
3091  for (i = 0; i <= fNpx; i++) {
3092  xbins[i] = gPad->PadtoX(xlogmin + i * dlogx);
3093  }
3094  histogram = new TH1D("Func", GetTitle(), fNpx, xbins);
3095  histogram->SetBit(TH1::kLogX);
3096  delete [] xbins;
3097  } else {
3098  histogram = new TH1D("Func", GetTitle(), fNpx, xmin, xmax);
3099  }
3100  if (fMinimum != -1111) histogram->SetMinimum(fMinimum);
3101  if (fMaximum != -1111) histogram->SetMaximum(fMaximum);
3102  histogram->SetDirectory(0);
3103  }
3104  R__ASSERT(histogram);
3105 
3106  // Restore axis titles.
3107  histogram->GetXaxis()->SetTitle(xtitle.Data());
3108  histogram->GetYaxis()->SetTitle(ytitle.Data());
3109  Double_t *parameters = GetParameters();
3110 
3111  InitArgs(xv, parameters);
3112  for (i = 1; i <= fNpx; i++) {
3113  xv[0] = histogram->GetBinCenter(i);
3114  histogram->SetBinContent(i, EvalPar(xv, parameters));
3115  }
3116 
3117  // Copy Function attributes to histogram attributes.
3118  histogram->SetBit(TH1::kNoStats);
3119  histogram->SetLineColor(GetLineColor());
3120  histogram->SetLineStyle(GetLineStyle());
3121  histogram->SetLineWidth(GetLineWidth());
3122  histogram->SetFillColor(GetFillColor());
3123  histogram->SetFillStyle(GetFillStyle());
3124  histogram->SetMarkerColor(GetMarkerColor());
3125  histogram->SetMarkerStyle(GetMarkerStyle());
3126  histogram->SetMarkerSize(GetMarkerSize());
3127 
3128  // update saved histogram in case it was deleted or if it is the first time the method is called
3129  // for example when called from TF1::GetHistogram()
3130  if (!fHistogram) fHistogram = histogram;
3131  return histogram;
3132 
3133 }
3134 
3135 
3136 ////////////////////////////////////////////////////////////////////////////////
3137 /// Release parameter number ipar If used in a fit, the parameter
3138 /// can vary freely. The parameter limits are reset to 0,0.
3139 
3141 {
3142  if (ipar < 0 || ipar > GetNpar() - 1) return;
3143  SetParLimits(ipar, 0, 0);
3144 }
3145 
3146 
3147 ////////////////////////////////////////////////////////////////////////////////
3148 /// Save values of function in array fSave
3149 
3151 {
3152  Double_t *parameters = GetParameters();
3153  //if (fSave != 0) {delete [] fSave; fSave = 0;}
3154  if (fParent && fParent->InheritsFrom(TH1::Class())) {
3155  //if parent is a histogram save the function at the center of the bins
3156  if ((xmin > 0 && xmax > 0) && TMath::Abs(TMath::Log10(xmax / xmin) > TMath::Log10(fNpx))) {
3157  TH1 *h = (TH1 *)fParent;
3158  Int_t bin1 = h->GetXaxis()->FindBin(xmin);
3159  Int_t bin2 = h->GetXaxis()->FindBin(xmax);
3160  int fNsave = bin2 - bin1 + 4;
3161  //fSave = new Double_t[fNsave];
3162  fSave.resize(fNsave);
3163  Double_t xv[1];
3164 
3165  InitArgs(xv, parameters);
3166  for (Int_t i = bin1; i <= bin2; i++) {
3167  xv[0] = h->GetXaxis()->GetBinCenter(i);
3168  fSave[i - bin1] = EvalPar(xv, parameters);
3169  }
3170  fSave[fNsave - 3] = xmin;
3171  fSave[fNsave - 2] = xmax;
3172  fSave[fNsave - 1] = xmax;
3173  return;
3174  }
3175  }
3176  int fNsave = fNpx + 3;
3177  if (fNsave <= 3) {
3178  return;
3179  }
3180  //fSave = new Double_t[fNsave];
3181  fSave.resize(fNsave);
3182  Double_t dx = (xmax - xmin) / fNpx;
3183  if (dx <= 0) {
3184  dx = (fXmax - fXmin) / fNpx;
3185  fNsave--;
3186  xmin = fXmin + 0.5 * dx;
3187  xmax = fXmax - 0.5 * dx;
3188  }
3189  Double_t xv[1];
3190  InitArgs(xv, parameters);
3191  for (Int_t i = 0; i <= fNpx; i++) {
3192  xv[0] = xmin + dx * i;
3193  fSave[i] = EvalPar(xv, parameters);
3194  }
3195  fSave[fNpx + 1] = xmin;
3196  fSave[fNpx + 2] = xmax;
3197 }
3198 
3199 
3200 ////////////////////////////////////////////////////////////////////////////////
3201 /// Save primitive as a C++ statement(s) on output stream out
3202 
3203 void TF1::SavePrimitive(std::ostream &out, Option_t *option /*= ""*/)
3204 {
3205  Int_t i;
3206  char quote = '"';
3207 
3208  // Save the function as C code independant from ROOT.
3209  if (strstr(option, "cc")) {
3210  out << "double " << GetName() << "(double xv) {" << std::endl;
3211  Double_t dx = (fXmax - fXmin) / (fNpx - 1);
3212  out << " double x[" << fNpx << "] = {" << std::endl;
3213  out << " ";
3214  Int_t n = 0;
3215  for (i = 0; i < fNpx; i++) {
3216  out << fXmin + dx *i ;
3217  if (i < fNpx - 1) out << ", ";
3218  if (n++ == 10) {
3219  out << std::endl;
3220  out << " ";
3221  n = 0;
3222  }
3223  }
3224  out << std::endl;
3225  out << " };" << std::endl;
3226  out << " double y[" << fNpx << "] = {" << std::endl;
3227  out << " ";
3228  n = 0;
3229  for (i = 0; i < fNpx; i++) {
3230  out << Eval(fXmin + dx * i);
3231  if (i < fNpx - 1) out << ", ";
3232  if (n++ == 10) {
3233  out << std::endl;
3234  out << " ";
3235  n = 0;
3236  }
3237  }
3238  out << std::endl;
3239  out << " };" << std::endl;
3240  out << " if (xv<x[0]) return y[0];" << std::endl;
3241  out << " if (xv>x[" << fNpx - 1 << "]) return y[" << fNpx - 1 << "];" << std::endl;
3242  out << " int i, j=0;" << std::endl;
3243  out << " for (i=1; i<" << fNpx << "; i++) { if (xv < x[i]) break; j++; }" << std::endl;
3244  out << " return y[j] + (y[j + 1] - y[j]) / (x[j + 1] - x[j]) * (xv - x[j]);" << std::endl;
3245  out << "}" << std::endl;
3246  return;
3247  }
3248 
3249  out << " " << std::endl;
3250 
3251  // Either f1Number is computed locally or set from outside
3252  static Int_t f1Number = 0;
3253  TString f1Name(GetName());
3254  const char *l = strstr(option, "#");
3255  if (l != 0) {
3256  sscanf(&l[1], "%d", &f1Number);
3257  } else {
3258  ++f1Number;
3259  }
3260  f1Name += f1Number;
3261 
3262  const char *addToGlobList = fParent ? ", TF1::EAddToList::kNo" : ", TF1::EAddToList::kDefault";
3263 
3264  if (!fType) {
3265  out << " TF1 *" << f1Name.Data() << " = new TF1(" << quote << GetName() << quote << "," << quote << GetTitle() << quote << "," << fXmin << "," << fXmax << addToGlobList << ");" << std::endl;
3266  if (fNpx != 100) {
3267  out << " " << f1Name.Data() << "->SetNpx(" << fNpx << ");" << std::endl;
3268  }
3269  } else {
3270  out << " TF1 *" << f1Name.Data() << " = new TF1(" << quote << "*" << GetName() << quote << "," << fXmin << "," << fXmax << "," << GetNpar() << ");" << std::endl;
3271  out << " //The original function : " << GetTitle() << " had originally been created by:" << std::endl;
3272  out << " //TF1 *" << GetName() << " = new TF1(" << quote << GetName() << quote << "," << GetTitle() << "," << fXmin << "," << fXmax << "," << GetNpar();
3273  out << ", 1" << addToGlobList << ");" << std::endl;
3274  out << " " << f1Name.Data() << "->SetRange(" << fXmin << "," << fXmax << ");" << std::endl;
3275  out << " " << f1Name.Data() << "->SetName(" << quote << GetName() << quote << ");" << std::endl;
3276  out << " " << f1Name.Data() << "->SetTitle(" << quote << GetTitle() << quote << ");" << std::endl;
3277  if (fNpx != 100) {
3278  out << " " << f1Name.Data() << "->SetNpx(" << fNpx << ");" << std::endl;
3279  }
3280  Double_t dx = (fXmax - fXmin) / fNpx;
3281  Double_t xv[1];
3282  Double_t *parameters = GetParameters();
3283  InitArgs(xv, parameters);
3284  for (i = 0; i <= fNpx; i++) {
3285  xv[0] = fXmin + dx * i;
3286  Double_t save = EvalPar(xv, parameters);
3287  out << " " << f1Name.Data() << "->SetSavedPoint(" << i << "," << save << ");" << std::endl;
3288  }
3289  out << " " << f1Name.Data() << "->SetSavedPoint(" << fNpx + 1 << "," << fXmin << ");" << std::endl;
3290  out << " " << f1Name.Data() << "->SetSavedPoint(" << fNpx + 2 << "," << fXmax << ");" << std::endl;
3291  }
3292 
3293  if (TestBit(kNotDraw)) {
3294  out << " " << f1Name.Data() << "->SetBit(TF1::kNotDraw);" << std::endl;
3295  }
3296  if (GetFillColor() != 0) {
3297  if (GetFillColor() > 228) {
3299  out << " " << f1Name.Data() << "->SetFillColor(ci);" << std::endl;
3300  } else
3301  out << " " << f1Name.Data() << "->SetFillColor(" << GetFillColor() << ");" << std::endl;
3302  }
3303  if (GetFillStyle() != 1001) {
3304  out << " " << f1Name.Data() << "->SetFillStyle(" << GetFillStyle() << ");" << std::endl;
3305  }
3306  if (GetMarkerColor() != 1) {
3307  if (GetMarkerColor() > 228) {
3309  out << " " << f1Name.Data() << "->SetMarkerColor(ci);" << std::endl;
3310  } else
3311  out << " " << f1Name.Data() << "->SetMarkerColor(" << GetMarkerColor() << ");" << std::endl;
3312  }
3313  if (GetMarkerStyle() != 1) {
3314  out << " " << f1Name.Data() << "->SetMarkerStyle(" << GetMarkerStyle() << ");" << std::endl;
3315  }
3316  if (GetMarkerSize() != 1) {
3317  out << " " << f1Name.Data() << "->SetMarkerSize(" << GetMarkerSize() << ");" << std::endl;
3318  }
3319  if (GetLineColor() != 1) {
3320  if (GetLineColor() > 228) {
3322  out << " " << f1Name.Data() << "->SetLineColor(ci);" << std::endl;
3323  } else
3324  out << " " << f1Name.Data() << "->SetLineColor(" << GetLineColor() << ");" << std::endl;
3325  }
3326  if (GetLineWidth() != 4) {
3327  out << " " << f1Name.Data() << "->SetLineWidth(" << GetLineWidth() << ");" << std::endl;
3328  }
3329  if (GetLineStyle() != 1) {
3330  out << " " << f1Name.Data() << "->SetLineStyle(" << GetLineStyle() << ");" << std::endl;
3331  }
3332  if (GetChisquare() != 0) {
3333  out << " " << f1Name.Data() << "->SetChisquare(" << GetChisquare() << ");" << std::endl;
3334  out << " " << f1Name.Data() << "->SetNDF(" << GetNDF() << ");" << std::endl;
3335  }
3336 
3337  if (GetXaxis()) GetXaxis()->SaveAttributes(out, f1Name.Data(), "->GetXaxis()");
3338  if (GetYaxis()) GetYaxis()->SaveAttributes(out, f1Name.Data(), "->GetYaxis()");
3339 
3340  Double_t parmin, parmax;
3341  for (i = 0; i < GetNpar(); i++) {
3342  out << " " << f1Name.Data() << "->SetParameter(" << i << "," << GetParameter(i) << ");" << std::endl;
3343  out << " " << f1Name.Data() << "->SetParError(" << i << "," << GetParError(i) << ");" << std::endl;
3344  GetParLimits(i, parmin, parmax);
3345  out << " " << f1Name.Data() << "->SetParLimits(" << i << "," << parmin << "," << parmax << ");" << std::endl;
3346  }
3347  if (!strstr(option, "nodraw")) {
3348  out << " " << f1Name.Data() << "->Draw("
3349  << quote << option << quote << ");" << std::endl;
3350  }
3351 }
3352 
3353 
3354 ////////////////////////////////////////////////////////////////////////////////
3355 /// Static function setting the current function.
3356 /// the current function may be accessed in static C-like functions
3357 /// when fitting or painting a function.
3358 
3360 {
3361  fgCurrent = f1;
3362 }
3363 
3364 ////////////////////////////////////////////////////////////////////////////////
3365 /// Set the result from the fit
3366 /// parameter values, errors, chi2, etc...
3367 /// Optionally a pointer to a vector (with size fNpar) of the parameter indices in the FitResult can be passed
3368 /// This is useful in the case of a combined fit with different functions, and the FitResult contains the global result
3369 /// By default it is assume that indpar = {0,1,2,....,fNpar-1}.
3370 
3371 void TF1::SetFitResult(const ROOT::Fit::FitResult &result, const Int_t *indpar)
3372 {
3373  Int_t npar = GetNpar();
3374  if (result.IsEmpty()) {
3375  Warning("SetFitResult", "Empty Fit result - nothing is set in TF1");
3376  return;
3377  }
3378  if (indpar == 0 && npar != (int) result.NPar()) {
3379  Error("SetFitResult", "Invalid Fit result passed - number of parameter is %d , different than TF1::GetNpar() = %d", npar, result.NPar());
3380  return;
3381  }
3382  if (result.Chi2() > 0)
3383  SetChisquare(result.Chi2());
3384  else
3385  SetChisquare(result.MinFcnValue());
3386 
3387  SetNDF(result.Ndf());
3388  SetNumberFitPoints(result.Ndf() + result.NFreeParameters());
3389 
3390 
3391  for (Int_t i = 0; i < npar; ++i) {
3392  Int_t ipar = (indpar != 0) ? indpar[i] : i;
3393  if (ipar < 0) continue;
3394  GetParameters()[i] = result.Parameter(ipar);
3395  // in case errors are not present do not set them
3396  if (ipar < (int) result.Errors().size())
3397  fParErrors[i] = result.Error(ipar);
3398  }
3399  //invalidate cached integral since parameters have changed
3400  Update();
3401 
3402 }
3403 
3404 
3405 ////////////////////////////////////////////////////////////////////////////////
3406 /// Set the maximum value along Y for this function
3407 /// In case the function is already drawn, set also the maximum in the
3408 /// helper histogram
3409 
3411 {
3412  fMaximum = maximum;
3413  if (fHistogram) fHistogram->SetMaximum(maximum);
3414  if (gPad) gPad->Modified();
3415 }
3416 
3417 
3418 ////////////////////////////////////////////////////////////////////////////////
3419 /// Set the minimum value along Y for this function
3420 /// In case the function is already drawn, set also the minimum in the
3421 /// helper histogram
3422 
3424 {
3425  fMinimum = minimum;
3426  if (fHistogram) fHistogram->SetMinimum(minimum);
3427  if (gPad) gPad->Modified();
3428 }
3429 
3430 
3431 ////////////////////////////////////////////////////////////////////////////////
3432 /// Set the number of degrees of freedom
3433 /// ndf should be the number of points used in a fit - the number of free parameters
3434 
3436 {
3437  fNDF = ndf;
3438 }
3439 
3440 
3441 ////////////////////////////////////////////////////////////////////////////////
3442 /// Set the number of points used to draw the function
3443 ///
3444 /// The default number of points along x is 100 for 1-d functions and 30 for 2-d/3-d functions
3445 /// You can increase this value to get a better resolution when drawing
3446 /// pictures with sharp peaks or to get a better result when using TF1::GetRandom
3447 /// the minimum number of points is 4, the maximum is 10000000 for 1-d and 10000 for 2-d/3-d functions
3448 
3450 {
3451  const Int_t minPx = 4;
3452  Int_t maxPx = 10000000;
3453  if (GetNdim() > 1) maxPx = 10000;
3454  if (npx >= minPx && npx <= maxPx) {
3455  fNpx = npx;
3456  } else {
3457  if (npx < minPx) fNpx = minPx;
3458  if (npx > maxPx) fNpx = maxPx;
3459  Warning("SetNpx", "Number of points must be >=%d && <= %d, fNpx set to %d", minPx, maxPx, fNpx);
3460  }
3461  Update();
3462 }
3463 ////////////////////////////////////////////////////////////////////////////////
3464 /// Set name of parameter number ipar
3465 
3466 void TF1::SetParName(Int_t ipar, const char *name)
3467 {
3468  if (fFormula) {
3469  if (ipar < 0 || ipar >= GetNpar()) return;
3470  fFormula->SetParName(ipar, name);
3471  } else
3472  fParams->SetParName(ipar, name);
3473 }
3474 
3475 ////////////////////////////////////////////////////////////////////////////////
3476 /// Set up to 10 parameter names
3477 
3478 void TF1::SetParNames(const char *name0, const char *name1, const char *name2, const char *name3, const char *name4,
3479  const char *name5, const char *name6, const char *name7, const char *name8, const char *name9, const char *name10)
3480 {
3481  if (fFormula)
3482  fFormula->SetParNames(name0, name1, name2, name3, name4, name5, name6, name7, name8, name9, name10);
3483  else
3484  fParams->SetParNames(name0, name1, name2, name3, name4, name5, name6, name7, name8, name9, name10);
3485 }
3486 ////////////////////////////////////////////////////////////////////////////////
3487 /// Set error for parameter number ipar
3488 
3490 {
3491  if (ipar < 0 || ipar > GetNpar() - 1) return;
3492  fParErrors[ipar] = error;
3493 }
3494 
3495 
3496 ////////////////////////////////////////////////////////////////////////////////
3497 /// Set errors for all active parameters
3498 /// when calling this function, the array errors must have at least fNpar values
3499 
3500 void TF1::SetParErrors(const Double_t *errors)
3501 {
3502  if (!errors) return;
3503  for (Int_t i = 0; i < GetNpar(); i++) fParErrors[i] = errors[i];
3504 }
3505 
3506 
3507 ////////////////////////////////////////////////////////////////////////////////
3508 /// Set limits for parameter ipar.
3509 ///
3510 /// The specified limits will be used in a fit operation
3511 /// when the option "B" is specified (Bounds).
3512 /// To fix a parameter, use TF1::FixParameter
3513 
3514 void TF1::SetParLimits(Int_t ipar, Double_t parmin, Double_t parmax)
3515 {
3516  Int_t npar = GetNpar();
3517  if (ipar < 0 || ipar > npar - 1) return;
3518  if (int(fParMin.size()) != npar) {
3519  fParMin.resize(npar);
3520  }
3521  if (int(fParMax.size()) != npar) {
3522  fParMax.resize(npar);
3523  }
3524  fParMin[ipar] = parmin;
3525  fParMax[ipar] = parmax;
3526 }
3527 
3528 
3529 ////////////////////////////////////////////////////////////////////////////////
3530 /// Initialize the upper and lower bounds to draw the function.
3531 ///
3532 /// The function range is also used in an histogram fit operation
3533 /// when the option "R" is specified.
3534 
3536 {
3537  fXmin = xmin;
3538  fXmax = xmax;
3539  if (fType == EFType::kCompositionFcn && fComposition) {
3540  fComposition->SetRange(xmin, xmax); // automatically updates sub-functions
3541  }
3542  Update();
3543 }
3544 
3545 
3546 ////////////////////////////////////////////////////////////////////////////////
3547 /// Restore value of function saved at point
3548 
3550 {
3551  if (fSave.size() == 0) {
3552  fSave.resize(fNpx + 3);
3553  }
3554  if (point < 0 || point >= int(fSave.size())) return;
3555  fSave[point] = value;
3556 }
3557 
3558 
3559 ////////////////////////////////////////////////////////////////////////////////
3560 /// Set function title
3561 /// if title has the form "fffffff;xxxx;yyyy", it is assumed that
3562 /// the function title is "fffffff" and "xxxx" and "yyyy" are the
3563 /// titles for the X and Y axis respectively.
3564 
3565 void TF1::SetTitle(const char *title)
3566 {
3567  if (!title) return;
3568  fTitle = title;
3569  if (!fHistogram) return;
3570  fHistogram->SetTitle(title);
3571  if (gPad) gPad->Modified();
3572 }
3573 
3574 
3575 ////////////////////////////////////////////////////////////////////////////////
3576 /// Stream a class object.
3577 
3578 void TF1::Streamer(TBuffer &b)
3579 {
3580  if (b.IsReading()) {
3581  UInt_t R__s, R__c;
3582  Version_t v = b.ReadVersion(&R__s, &R__c);
3583  // process new version with new TFormula class which is contained in TF1
3584  //printf("reading TF1....- version %d..\n",v);
3585 
3586  if (v > 7) {
3587  // new classes with new TFormula
3588  // need to register the objects
3589  b.ReadClassBuffer(TF1::Class(), this, v, R__s, R__c);
3590  if (!TestBit(kNotGlobal)) {
3592  gROOT->GetListOfFunctions()->Add(this);
3593  }
3594  if (v >= 10)
3595  fComposition = std::unique_ptr<TF1AbsComposition>(fComposition_ptr);
3596  return;
3597  } else {
3598  ROOT::v5::TF1Data fold;
3599  //printf("Reading TF1 as v5::TF1Data- version %d \n",v);
3600  fold.Streamer(b, v, R__s, R__c, TF1::Class());
3601  // convert old TF1 to new one
3602  ((TF1v5Convert *)this)->Convert(fold);
3603  }
3604  }
3605 
3606  // Writing
3607  else {
3608  Int_t saved = 0;
3609  // save not-formula functions as array of points
3610  if (fType > 0 && fSave.empty() && fType != EFType::kCompositionFcn) {
3611  saved = 1;
3612  Save(fXmin, fXmax, 0, 0, 0, 0);
3613  }
3614  if (fType == EFType::kCompositionFcn)
3616  else
3617  fComposition_ptr = nullptr;
3618  b.WriteClassBuffer(TF1::Class(), this);
3619 
3620  // clear vector contents
3621  if (saved) {
3622  fSave.clear();
3623  }
3624  }
3625 }
3626 
3627 
3628 ////////////////////////////////////////////////////////////////////////////////
3629 /// Called by functions such as SetRange, SetNpx, SetParameters
3630 /// to force the deletion of the associated histogram or Integral
3631 
3633 {
3634  delete fHistogram;
3635  fHistogram = 0;
3636  if (!fIntegral.empty()) {
3637  fIntegral.clear();
3638  fAlpha.clear();
3639  fBeta.clear();
3640  fGamma.clear();
3641  }
3642  if (fNormalized) {
3643  // need to compute the integral of the not-normalized function
3644  fNormalized = false;
3645  fNormIntegral = Integral(fXmin, fXmax, 0.0);
3646  fNormalized = true;
3647  } else
3648  fNormIntegral = 0;
3649 
3650  // std::vector<double>x(fNdim);
3651  // if ((fType == 1) && !fFunctor->Empty()) (*fFunctor)x.data(), (Double_t*)fParams);
3652  if (fType == EFType::kCompositionFcn && fComposition) {
3653  // double-check that the parameters are correct
3654  fComposition->SetParameters(GetParameters());
3655 
3656  fComposition->Update(); // should not be necessary, but just to be safe
3657  }
3658 }
3659 
3660 ////////////////////////////////////////////////////////////////////////////////
3661 /// Static function to set the global flag to reject points
3662 /// the fgRejectPoint global flag is tested by all fit functions
3663 /// if TRUE the point is not included in the fit.
3664 /// This flag can be set by a user in a fitting function.
3665 /// The fgRejectPoint flag is reset by the TH1 and TGraph fitting functions.
3666 
3668 {
3669  fgRejectPoint = reject;
3670 }
3671 
3672 
3673 ////////////////////////////////////////////////////////////////////////////////
3674 /// See TF1::RejectPoint above
3675 
3677 {
3678  return fgRejectPoint;
3679 }
3680 
3681 ////////////////////////////////////////////////////////////////////////////////
3682 /// Return nth moment of function between a and b
3683 ///
3684 /// See TF1::Integral() for parameter definitions
3685 
3687 {
3688  // wrapped function in interface for integral calculation
3689  // using abs value of integral
3690 
3691  TF1_EvalWrapper func(this, params, kTRUE, n);
3692 
3694 
3695  giod.SetFunction(func);
3696  giod.SetRelTolerance(epsilon);
3697 
3698  Double_t norm = giod.Integral(a, b);
3699  if (norm == 0) {
3700  Error("Moment", "Integral zero over range");
3701  return 0;
3702  }
3703 
3704  // calculate now integral of x^n f(x)
3705  // wrapped the member function EvalNum in interface required by integrator using the functor class
3706  ROOT::Math::Functor1D xnfunc(&func, &TF1_EvalWrapper::EvalNMom);
3707  giod.SetFunction(xnfunc);
3708 
3709  Double_t res = giod.Integral(a, b) / norm;
3710 
3711  return res;
3712 }
3713 
3714 
3715 ////////////////////////////////////////////////////////////////////////////////
3716 /// Return nth central moment of function between a and b
3717 /// (i.e the n-th moment around the mean value)
3718 ///
3719 /// See TF1::Integral() for parameter definitions
3720 ///
3721 /// \author Gene Van Buren <gene@bnl.gov>
3722 
3724 {
3725  TF1_EvalWrapper func(this, params, kTRUE, n);
3726 
3728 
3729  giod.SetFunction(func);
3730  giod.SetRelTolerance(epsilon);
3731 
3732  Double_t norm = giod.Integral(a, b);
3733  if (norm == 0) {
3734  Error("Moment", "Integral zero over range");
3735  return 0;
3736  }
3737 
3738  // calculate now integral of xf(x)
3739  // wrapped the member function EvalFirstMom in interface required by integrator using the functor class
3740  ROOT::Math::Functor1D xfunc(&func, &TF1_EvalWrapper::EvalFirstMom);
3741  giod.SetFunction(xfunc);
3742 
3743  // estimate of mean value
3744  Double_t xbar = giod.Integral(a, b) / norm;
3745 
3746  // use different mean value in function wrapper
3747  func.fX0 = xbar;
3748  ROOT::Math::Functor1D xnfunc(&func, &TF1_EvalWrapper::EvalNMom);
3749  giod.SetFunction(xnfunc);
3750 
3751  Double_t res = giod.Integral(a, b) / norm;
3752  return res;
3753 }
3754 
3755 
3756 //______________________________________________________________________________
3757 // some useful static utility functions to compute sampling points for IntegralFast
3758 ////////////////////////////////////////////////////////////////////////////////
3759 /// Type safe interface (static method)
3760 /// The number of sampling points are taken from the TGraph
3761 
3762 #ifdef INTHEFUTURE
3764 {
3765  if (!g) return;
3766  CalcGaussLegendreSamplingPoints(g->GetN(), g->GetX(), g->GetY(), eps);
3767 }
3768 
3769 
3770 ////////////////////////////////////////////////////////////////////////////////
3771 /// Type safe interface (static method)
3772 /// A TGraph is created with new with num points and the pointer to the
3773 /// graph is returned by the function. It is the responsibility of the
3774 /// user to delete the object.
3775 /// if num is invalid (<=0) NULL is returned
3776 
3778 {
3779  if (num <= 0)
3780  return 0;
3781 
3782  TGraph *g = new TGraph(num);
3783  CalcGaussLegendreSamplingPoints(g->GetN(), g->GetX(), g->GetY(), eps);
3784  return g;
3785 }
3786 #endif
3787 
3788 
3789 ////////////////////////////////////////////////////////////////////////////////
3790 /// Type: unsafe but fast interface filling the arrays x and w (static method)
3791 ///
3792 /// Given the number of sampling points this routine fills the arrays x and w
3793 /// of length num, containing the abscissa and weight of the Gauss-Legendre
3794 /// n-point quadrature formula.
3795 ///
3796 /// Gauss-Legendre:
3797 /** \f[
3798  W(x)=1 -1<x<1 \\
3799  (j+1)P_{j+1} = (2j+1)xP_j-jP_{j-1}
3800  \f]
3801 **/
3802 /// num is the number of sampling points (>0)
3803 /// x and w are arrays of size num
3804 /// eps is the relative precision
3805 ///
3806 /// If num<=0 or eps<=0 no action is done.
3807 ///
3808 /// Reference: Numerical Recipes in C, Second Edition
3809 
3811 {
3812  // This function is just kept like this for backward compatibility!
3813 
3815  gli.GetWeightVectors(x, w);
3816 
3817 
3818 }
3819 
3820 
3821 /** \class TF1Parameters
3822 TF1 Parameters class
3823 */
3824 
3825 ////////////////////////////////////////////////////////////////////////////////
3826 /// Returns the parameter number given a name
3827 /// not very efficient but list of parameters is typically small
3828 /// could use a map if needed
3829 
3831 {
3832  for (unsigned int i = 0; i < fParNames.size(); ++i) {
3833  if (fParNames[i] == std::string(name)) return i;
3834  }
3835  return -1;
3836 }
3837 
3838 ////////////////////////////////////////////////////////////////////////////////
3839 /// Set parameter values
3840 
3842  Double_t p5, Double_t p6, Double_t p7, Double_t p8,
3843  Double_t p9, Double_t p10)
3844 {
3845  unsigned int npar = fParameters.size();
3846  if (npar > 0) fParameters[0] = p0;
3847  if (npar > 1) fParameters[1] = p1;
3848  if (npar > 2) fParameters[2] = p2;
3849  if (npar > 3) fParameters[3] = p3;
3850  if (npar > 4) fParameters[4] = p4;
3851  if (npar > 5) fParameters[5] = p5;
3852  if (npar > 6) fParameters[6] = p6;
3853  if (npar > 7) fParameters[7] = p7;
3854  if (npar > 8) fParameters[8] = p8;
3855  if (npar > 9) fParameters[9] = p9;
3856  if (npar > 10) fParameters[10] = p10;
3857 }
3858 
3859 ////////////////////////////////////////////////////////////////////////////////
3860 /// Set parameter names
3861 
3862 void TF1Parameters::SetParNames(const char *name0, const char *name1, const char *name2, const char *name3,
3863  const char *name4, const char *name5, const char *name6, const char *name7,
3864  const char *name8, const char *name9, const char *name10)
3865 {
3866  unsigned int npar = fParNames.size();
3867  if (npar > 0) fParNames[0] = name0;
3868  if (npar > 1) fParNames[1] = name1;
3869  if (npar > 2) fParNames[2] = name2;
3870  if (npar > 3) fParNames[3] = name3;
3871  if (npar > 4) fParNames[4] = name4;
3872  if (npar > 5) fParNames[5] = name5;
3873  if (npar > 6) fParNames[6] = name6;
3874  if (npar > 7) fParNames[7] = name7;
3875  if (npar > 8) fParNames[8] = name8;
3876  if (npar > 9) fParNames[9] = name9;
3877  if (npar > 10) fParNames[10] = name10;
3878 }
TString fTitle
Definition: TNamed.h:33
virtual const char * GetName() const
Returns name of object.
Definition: TNamed.h:47
Double_t fMaximum
Definition: TF1.h:250
virtual void Print(Option_t *option="") const
Print some global quantities for this histogram.
Definition: TH1.cxx:6620
virtual Double_t GetMaximumStored() const
Definition: TH1.h:284
virtual Int_t GetQuantiles(Int_t nprobSum, Double_t *q, const Double_t *probSum)
Compute Quantiles for density distribution of this function.
Definition: TF1.cxx:1990
virtual void SetLineWidth(Width_t lwidth)
Set the line width.
Definition: TAttLine.h:43
virtual Double_t GetMaximum(Double_t xmin=0, Double_t xmax=0, Double_t epsilon=1.E-10, Int_t maxiter=100, Bool_t logx=false) const
Returns the maximum value of the function.
Definition: TF1.cxx:1611
double Derivative3(double x)
Returns the third derivative of the function at point x, computed by Richardson&#39;s extrapolation metho...
virtual Double_t GetMaximum(Double_t maxval=FLT_MAX) const
Return maximum value smaller than maxval of bins in the range, unless the value has been overridden b...
Definition: TH1.cxx:7984
virtual void SaveAttributes(std::ostream &out, const char *name, const char *subname)
Save axis attributes as C++ statement(s) on output stream out.
Definition: TAxis.cxx:658
virtual void Paint(Option_t *option="")
Control routine to paint any kind of histograms.
Definition: TH1.cxx:5815
An array of TObjects.
Definition: TObjArray.h:37
virtual void SetParameters(const Double_t *params)
Definition: TF1.h:638
virtual Double_t GetRandom()
Return a random number following this function shape.
Definition: TF1.cxx:2091
Int_t fNpx
Definition: TF1.h:244
Bool_t fNormalized
Pointer to MethodCall in case of interpreted function.
Definition: TF1.h:262
double IntegralUp(double a)
Returns Integral of function on an upper semi-infinite interval.
float xmin
Definition: THbookFile.cxx:93
virtual Double_t GetBinCenter(Int_t bin) const
Return bin center for 1D histogram.
Definition: TH1.cxx:8575
virtual void Info(const char *method, const char *msgfmt,...) const
Issue info message.
Definition: TObject.cxx:861
Double_t Eval(Double_t x) const
Param Functor class for Multidimensional functions.
Definition: ParamFunctor.h:274
double Error() const
Return the estimate of the absolute Error of the last Integral calculation.
Interface (abstract class) for generic functions objects of one-dimension Provides a method to evalua...
Definition: IFunction.h:135
virtual void SetNpx(Int_t npx=100)
Set the number of points used to draw the function.
Definition: TF1.cxx:3449
virtual void ReleaseParameter(Int_t ipar)
Release parameter number ipar If used in a fit, the parameter can vary freely.
Definition: TF1.cxx:3140
long long Long64_t
Definition: RtypesCore.h:69
virtual void SetMaximum(Double_t maximum=-1111)
Definition: TH1.h:394
bool Solve(int maxIter=100, double absTol=1E-8, double relTol=1E-10)
Returns the X value corresponding to the function value fy for (xmin<x<xmax).
virtual Double_t GetMinimumX(Double_t xmin=0, Double_t xmax=0, Double_t epsilon=1.E-10, Int_t maxiter=100, Bool_t logx=false) const
Returns the X value corresponding to the minimum value of the function on the (xmin, xmax) interval.
Definition: TF1.cxx:1820
auto * m
Definition: textangle.C:8
virtual void SetParName(Int_t ipar, const char *name)
Set name of parameter number ipar.
Definition: TF1.cxx:3466
virtual void Copy(TObject &f1) const
Copy this F1 to a new F1.
Definition: TF1.cxx:994
virtual void SetLimits(Double_t xmin, Double_t xmax)
Definition: TAxis.h:154
const char * GetParName(Int_t ipar) const
Returns the available number of logical cores.
Definition: StringConv.hxx:21
int Status() const
return the Error Status of the last Integral calculation
Definition: Integrator.h:417
short Version_t
Definition: RtypesCore.h:61
double Error(unsigned int i) const
parameter error by index
Definition: FitResult.h:186
virtual void SetParNames(const char *name0="p0", const char *name1="p1", const char *name2="p2", const char *name3="p3", const char *name4="p4", const char *name5="p5", const char *name6="p6", const char *name7="p7", const char *name8="p8", const char *name9="p9", const char *name10="p10")
Set up to 10 parameter names.
Definition: TF1.cxx:3478
Collectable string class.
Definition: TObjString.h:28
virtual void SetDirectory(TDirectory *dir)
By default when an histogram is created, it is added to the list of histogram objects in the current ...
Definition: TH1.cxx:8371
virtual TF1 * DrawCopy(Option_t *option="") const
Draw a copy of this function with its current attributes.
Definition: TF1.cxx:1364
static std::atomic< Bool_t > fgAbsValue
Definition: TF1.h:300
const char Option_t
Definition: RtypesCore.h:62
static void R__v5TF1Updater(Int_t nobjects, TObject **from, TObject **to)
Definition: TF1.cxx:133
float ymin
Definition: THbookFile.cxx:93
virtual Int_t GetNumberFreeParameters() const
Return the number of free parameters.
Definition: TF1.cxx:1897
const std::vector< double > & Errors() const
parameter errors (return st::vector)
Definition: FitResult.h:169
Double_t QuietNaN()
Returns a quiet NaN as defined by IEEE 754
Definition: TMath.h:891
#define g(i)
Definition: RSha256.hxx:105
EFType
Definition: TF1.h:231
unsigned int NPar() const
total number of parameters (abbreviation)
Definition: FitResult.h:131
virtual Double_t GetMinimumStored() const
Definition: TH1.h:288
TString & ReplaceAll(const TString &s1, const TString &s2)
Definition: TString.h:687
TObject * Clone(const char *newname=0) const
Make a complete copy of the underlying object.
Definition: TF1.cxx:1065
R__EXTERN TStyle * gStyle
Definition: TStyle.h:410
const Double_t * GetArray() const
Definition: TArrayD.h:43
Double_t fXmax
Definition: TF1Data.h:40
int Status() const
return status of integration
virtual void SetOwner(Bool_t enable=kTRUE)
Set whether this collection is the owner (enable==true) of its content.
void DoInitialize(EAddToList addToGlobList)
Common initialization of the TF1.
Definition: TF1.cxx:786
#define BIT(n)
Definition: Rtypes.h:83
virtual Double_t IntegralFast(Int_t num, Double_t *x, Double_t *w, Double_t a, Double_t b, Double_t *params=0, Double_t epsilon=1e-12)
Gauss-Legendre integral, see CalcGaussLegendreSamplingPoints.
Definition: TF1.cxx:2765
static TF1 * GetCurrent()
Static function returning the current function being processed.
Definition: TF1.cxx:1571
static void SaveColor(std::ostream &out, Int_t ci)
Save a color with index > 228 as a C++ statement(s) on output stream out.
Definition: TColor.cxx:2108
Int_t fNpar
Definition: TF1.h:242
virtual void Copy(TObject &f1) const
Copy this to obj.
Definition: TFormula.cxx:637
void SetParNames(const char *name0="p0", const char *name1="p1", const char *name2="p2", const char *name3="p3", const char *name4="p4", const char *name5="p5", const char *name6="p6", const char *name7="p7", const char *name8="p8", const char *name9="p9", const char *name10="p10")
TObject * fParent
Array gamma.
Definition: TF1.h:259
Double_t * fParErrors
Definition: TF1Data.h:47
virtual void SetRange(Double_t xmin, Double_t xmax)
Initialize the upper and lower bounds to draw the function.
Definition: TF1.cxx:3535
virtual void ExecuteEvent(Int_t event, Int_t px, Int_t py)
Execute action corresponding to one event.
Definition: TF1.cxx:1543
void ToUpper()
Change string to upper case.
Definition: TString.cxx:1138
Double_t fChisquare
Definition: TF1.h:248
Buffer base class used for serializing objects.
Definition: TBuffer.h:42
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:3150
#define R__ASSERT(e)
Definition: TError.h:96
virtual Double_t IntegralMultiple(Int_t n, const Double_t *a, const Double_t *b, Int_t maxpts, Double_t epsrel, Double_t epsabs, Double_t &relerr, Int_t &nfnevl, Int_t &ifail)
This function computes, to an attempted specified accuracy, the value of the integral.
Definition: TF1.cxx:2838
virtual void SetMinimum(Double_t minimum=-1111)
Definition: TH1.h:395
#define gROOT
Definition: TROOT.h:415
R__ALWAYS_INLINE Bool_t TestBit(UInt_t f) const
Definition: TObject.h:172
double Integral(const double *xmin, const double *xmax)
evaluate the integral with the previously given function between xmin[] and xmax[] ...
std::vector< std::string > fParNames
Definition: TF1.h:138
unsigned int Ndf() const
Number of degree of freedom.
Definition: FitResult.h:163
static ROOT::Math::Minimizer * CreateMinimizer(const std::string &minimizerType="", const std::string &algoType="")
static method to create the corrisponding Minimizer given the string Supported Minimizers types are: ...
Definition: Factory.cxx:63
Basic string class.
Definition: TString.h:131
Class to Wrap a ROOT Function class (like TF1) in a IParamFunction interface of one dimensions to be ...
Definition: WrappedTF1.h:37
Int_t GetNpar() const
Definition: TFormula.h:222
TAxis * GetXaxis() const
Get x axis of the function.
Definition: TF1.cxx:2399
virtual void SetNumberFitPoints(Int_t npfits)
Definition: TF1.h:618
static Bool_t DefaultAddToGlobalList(Bool_t on=kTRUE)
Static method to add/avoid to add automatically functions to the global list (gROOT->GetListOfFunctio...
Definition: TF1.cxx:822
static Double_t DerivativeError()
Static function returning the error of the last call to the of Derivative&#39;s functions.
Definition: TF1.cxx:1278
Double_t fNormIntegral
Definition: TF1.h:263
#define f(i)
Definition: RSha256.hxx:104
void ToLower()
Change string to lower-case.
Definition: TString.cxx:1125
int Int_t
Definition: RtypesCore.h:41
bool Bool_t
Definition: RtypesCore.h:59
virtual Double_t GetMaximumX(Double_t xmin=0, Double_t xmax=0, Double_t epsilon=1.E-10, Int_t maxiter=100, Bool_t logx=false) const
Returns the X value corresponding to the maximum value of the function.
Definition: TF1.cxx:1652
void(*)(Int_t nobjects, TObject **from, TObject **to) TF1Updater_t
Definition: TF1.cxx:58
R__EXTERN TVirtualMutex * gROOTMutex
Definition: TROOT.h:59
void Copy(TAttMarker &attmarker) const
Copy this marker attributes to a new TAttMarker.
Definition: TAttMarker.cxx:235
virtual void Draw(Option_t *option="")
Draw this function with its current attributes.
Definition: TF1.cxx:1334
virtual void SetFillStyle(Style_t fstyle)
Set the fill area style.
Definition: TAttFill.h:39
virtual Double_t GetParError(Int_t ipar) const
Return value of parameter number ipar.
Definition: TF1.cxx:1927
double Parameter(unsigned int i) const
parameter value by index
Definition: FitResult.h:181
virtual Double_t Derivative2(Double_t x, Double_t *params=0, Double_t epsilon=0.001) const
Returns the second derivative of the function at point x, computed by Richardson&#39;s extrapolation meth...
Definition: TF1.cxx:1179
Bool_t IsNaN(Double_t x)
Definition: TMath.h:882
double MinFcnValue() const
Return value of the objective function (chi2 or likelihood) used in the fit.
Definition: FitResult.h:120
double Integral(double a, double b)
Returns Integral of function between a and b.
virtual void SetMinimum(Double_t minimum=-1111)
Set the minimum value along Y for this function In case the function is already drawn, set also the minimum in the helper histogram.
Definition: TF1.cxx:3423
double Derivative2(double x)
Returns the second derivative of the function at point x, computed by Richardson&#39;s extrapolation meth...
TObject * At(Int_t idx) const
Definition: TObjArray.h:166
virtual Double_t Integral(Double_t a, Double_t b, Double_t epsrel=1.e-12)
IntegralOneDim or analytical integral.
Definition: TF1.cxx:2519
User class for performing function integration.
virtual TObject * DrawIntegral(Option_t *option="al")
Draw integral of this function.
Definition: TF1.cxx:1411
TRObject operator()(const T1 &t1) const
Short_t Abs(Short_t d)
Definition: TMathBase.h:120
virtual void Draw(Option_t *chopt="")
Draw this graph with its current attributes.
Definition: TGraph.cxx:751
bit set when zooming on Y axis
Definition: TH1.h:164
virtual TVirtualPad * cd(Int_t subpadnumber=0)=0
virtual Double_t Derivative3(Double_t x, Double_t *params=0, Double_t epsilon=0.001) const
Returns the third derivative of the function at point x, computed by Richardson&#39;s extrapolation metho...
Definition: TF1.cxx:1244
virtual Width_t GetLineWidth() const
Return the line width.
Definition: TAttLine.h:35
virtual Int_t GetNdim() const
Definition: TFormula.h:237
don&#39;t draw stats box
Definition: TH1.h:160
Double_t Prob(Double_t chi2, Int_t ndf)
Width_t GetFuncWidth() const
Definition: TStyle.h:210
LongDouble_t Power(LongDouble_t x, LongDouble_t y)
Definition: TMath.h:725
virtual bool SetLimitedVariable(unsigned int ivar, const std::string &name, double val, double step, double lower, double upper)
set a new upper/lower limited variable (override if minimizer supports them ) otherwise as default se...
Definition: Minimizer.h:163
void SetBit(UInt_t f, Bool_t set)
Set or unset the user status bits as specified in f.
Definition: TObject.cxx:693
double IntegralLow(double b)
Returns Integral of function on a lower semi-infinite interval.
double beta(double x, double y)
Calculates the beta function.
static void SetCurrent(TF1 *f1)
Static function setting the current function.
Definition: TF1.cxx:3359
if object in a list can be deleted
Definition: TObject.h:58
TF1 Parameters class.
Definition: TF1.h:47
virtual Double_t Derivative(Double_t x, Double_t *params=0, Double_t epsilon=0.001) const
Returns the first derivative of the function at point x, computed by Richardson&#39;s extrapolation metho...
Definition: TF1.cxx:1114
virtual void AppendPad(Option_t *option="")
Append graphics object to current pad.
Definition: TObject.cxx:104
virtual void RecursiveRemove(TObject *obj)
Recursively remove this object from a list.
Definition: TObject.cxx:571
TF1AbsComposition * fComposition_ptr
Pointer to composition (NSUM or CONV)
Definition: TF1.h:268
virtual Double_t GetProb() const
Return the fit probability.
Definition: TF1.cxx:1952
virtual Style_t GetMarkerStyle() const
Return the marker style.
Definition: TAttMarker.h:32
Template class to wrap any C++ callable object which takes one argument i.e.
void SetParamPtrs(void *paramArr, Int_t nparam=-1)
ParamArr is an array containing the function argument values.
Marker Attributes class.
Definition: TAttMarker.h:19
virtual Int_t DistancetoPrimitive(Int_t px, Int_t py)
Compute distance from point px,py to a function.
Definition: TF1.cxx:1294
virtual Style_t GetLineStyle() const
Return the line style.
Definition: TAttLine.h:34
virtual Int_t GetNDF() const
Return the number of degrees of freedom in the fit the fNDF parameter has been previously computed du...
Definition: TF1.cxx:1886
TAxis * GetYaxis() const
Get y axis of the function.
Definition: TF1.cxx:2410
virtual Double_t Moment(Double_t n, Double_t a, Double_t b, const Double_t *params=0, Double_t epsilon=0.000001)
Return nth moment of function between a and b.
Definition: TF1.cxx:3686
virtual void SetRelTolerance(double eps)
Set the desired relative Error.
EFType fType
Definition: TF1.h:245
Class wrapping convolution of two functions.
void SetParameters(const double *p)
set parameter values need to call also SetParameters in TF1 in ace some other operations (re-normaliz...
Definition: WrappedTF1.h:88
TFormula * fFormula
Functor object to wrap any C++ callable object.
Definition: TF1.h:265
Double_t GetXmin() const
Definition: TAxis.h:133
static const double x2[5]
Fill Area Attributes class.
Definition: TAttFill.h:19
Double_t x[n]
Definition: legend1.C:17
static TString Format(const char *fmt,...)
Static method which formats a string using a printf style format descriptor and return a TString...
Definition: TString.cxx:2311
void Class()
Definition: Class.C:29
Int_t fNpfits
Definition: TF1.h:246
Int_t fNdim
Definition: TF1.h:243
void Copy(TAttLine &attline) const
Copy this line attributes to a new TAttLine.
Definition: TAttLine.cxx:172
The TNamed class is the base class for all named ROOT classes.
Definition: TNamed.h:29
virtual void DrawF1(Double_t xmin, Double_t xmax, Option_t *option="")
Draw function between xmin and xmax.
Definition: TF1.cxx:1427
Abstract Minimizer class, defining the interface for the various minimizer (like Minuit2, Minuit, GSL, etc..) Plug-in&#39;s exist in ROOT to be able to instantiate the derived classes like ROOT::Math::GSLMinimizer or ROOT::Math::Minuit2Minimizer via the plug-in manager.
Definition: Minimizer.h:78
Double_t Log10(Double_t x)
Definition: TMath.h:754
virtual double MinValue() const =0
return minimum function value
int NEval() const
return number of function evaluations in calculating the integral
virtual Double_t GetBinCenter(Int_t bin) const
Return center of bin.
Definition: TAxis.cxx:475
TF1()
TF1 default constructor.
Definition: TF1.cxx:490
virtual bool Minimize()=0
method to perform the minimization
void SetParameters(const Double_t *params)
Definition: TF1.h:105
virtual void SetMarkerColor(Color_t mcolor=1)
Set the marker color.
Definition: TAttMarker.h:38
TString & Append(const char *cs)
Definition: TString.h:559
double Derivative1(double x)
Returns the first derivative of the function at point x, computed by Richardson&#39;s extrapolation metho...
virtual Size_t GetMarkerSize() const
Return the marker size.
Definition: TAttMarker.h:33
void function(const Char_t *name_, T fun, const Char_t *docstring=0)
Definition: RExports.h:151
static IntegrationOneDim::Type DefaultIntegratorType()
Ssiz_t First(char c) const
Find first occurrence of a character c.
Definition: TString.cxx:499
virtual Int_t GetNdim() const
Definition: TF1.h:479
you should not use this method at all Int_t Int_t Double_t bm
Definition: TRolke.cxx:630
virtual const double * X() const =0
return pointer to X values at the minimum
Method or function calling interface.
Definition: TMethodCall.h:37
User class for performing function minimization.
void SetParName(Int_t ipar, const char *name)
TVirtualPad is an abstract base class for the Pad and Canvas classes.
Definition: TVirtualPad.h:50
Double_t Infinity()
Returns an infinity as defined by the IEEE standard.
Definition: TMath.h:904
int TermCoeffLength(TString &term)
Definition: TF1.cxx:908
Int_t GetNdim() const
Definition: TFormula.h:221
double Error() const
return the estimate of the absolute Error of the last Integral calculation
Definition: Integrator.h:412
std::vector< Double_t > fIntegral
Definition: TF1.h:255
double Integral(const double *xmin, const double *xmax)
evaluate the integral with the previously given function between xmin[] and xmax[] ...
int Status() const
return the status of the last integration - 0 in case of success
virtual Double_t Rndm()
Machine independent random number generator.
Definition: TRandom.cxx:541
virtual Double_t GradientPar(Int_t ipar, const Double_t *x, Double_t eps=0.01)
Compute the gradient (derivative) wrt a parameter ipar.
Definition: TF1.cxx:2444
double IntegralLow(const IGenFunction &f, double b)
evaluate the Integral of a function f over the over the semi-infinite interval (-inf,b)
Definition: Integrator.h:292
virtual const char * GetParName(Int_t ipar) const
Definition: TF1.h:523
static TF1 * fgCurrent
Definition: TF1.h:303
Bool_t AreEqualRel(Double_t af, Double_t bf, Double_t relPrec)
Definition: TMath.h:418
virtual Int_t DistancetoPrimitive(Int_t px, Int_t py)
Compute distance from point px,py to a line.
Definition: TH1.cxx:2718
virtual void SetLineColor(Color_t lcolor)
Set the line color.
Definition: TAttLine.h:40
Double_t EvalPar(const Double_t *x, const Double_t *params=0) const
double Integral(Function &f, double a, double b)
evaluate the Integral of a function f over the defined interval (a,b)
Definition: Integrator.h:496
double Root() const
Returns root value.
Using a TBrowser one can browse all ROOT objects.
Definition: TBrowser.h:37
virtual TF1FunctorPointer * Clone() const =0
virtual void ExecuteEvent(Int_t event, Int_t px, Int_t py)
Execute action corresponding to one event.
Definition: TH1.cxx:3153
double gamma(double x)
virtual void SetChisquare(Double_t chi2)
Definition: TF1.h:606
float ymax
Definition: THbookFile.cxx:93
virtual void SetParLimits(Int_t ipar, Double_t parmin, Double_t parmax)
Set limits for parameter ipar.
Definition: TF1.cxx:3514
double IntegralUp(const IGenFunction &f, double a)
evaluate the Integral of a function f over the semi-infinite interval (a,+inf)
Definition: Integrator.h:274
TF1::EAddToList GetGlobalListOption(Option_t *opt)
Definition: TF1.cxx:672
static Bool_t fgRejectPoint
Definition: TF1.h:301
virtual void SetFunction(const ROOT::Math::IMultiGenFunction &func)=0
set the function to minimize
std::vector< Double_t > fParameters
Definition: TF1.h:137
TH1 * fHistogram
Parent object hooking this function (if one)
Definition: TF1.h:260
unsigned int NFreeParameters() const
get total number of free parameters
Definition: FitResult.h:134
virtual Double_t GetMinimum(Double_t xmin=0, Double_t xmax=0, Double_t epsilon=1.E-10, Int_t maxiter=100, Bool_t logx=false) const
Returns the minimum value of the function on the (xmin, xmax) interval.
Definition: TF1.cxx:1693
ROOT::R::TRInterface & r
Definition: Object.C:4
const char * GetTitle() const
Returns title of object.
Definition: TAxis.h:129
virtual Bool_t IsValid() const
Return kTRUE if the function is valid.
Definition: TF1.cxx:2869
Class to manage histogram axis.
Definition: TAxis.h:30
static const std::string & DefaultMinimizerType()
virtual void SetFillColor(Color_t fcolor)
Set the fill area color.
Definition: TAttFill.h:37
auto * a
Definition: textangle.C:12
void SetLogScan(bool on)
Set a log grid scan (default is equidistant bins) will work only if xlow > 0.
double IntegralError(TF1 *func, Int_t ndim, const double *a, const double *b, const double *params, const double *covmat, double epsilon)
Definition: TF1Helper.cxx:38
The Formula class.
Definition: TFormula.h:83
virtual void SetParError(Int_t ipar, Double_t error)
Set error for parameter number ipar.
Definition: TF1.cxx:3489
virtual char * GetObjectInfo(Int_t px, Int_t py) const
Redefines TObject::GetObjectInfo.
Definition: TF1.cxx:1915
virtual Bool_t InheritsFrom(const char *classname) const
Returns kTRUE if object inherits from class "classname".
Definition: TObject.cxx:442
virtual void SetBinContent(Int_t bin, Double_t content)
Set bin content see convention for numbering bins in TH1::GetBin In case the bin number is greater th...
Definition: TH1.cxx:8656
std::unique_ptr< TF1AbsComposition > fComposition
Definition: TF1.h:267
unsigned int UInt_t
Definition: RtypesCore.h:42
virtual void Error(const char *method, const char *msgfmt,...) const
Issue error message.
Definition: TObject.cxx:887
char * Form(const char *fmt,...)
User class for performing function integration.
virtual void FixParameter(Int_t ipar, Double_t value)
Fix the value of a parameter The specified value will be used in a fit operation. ...
Definition: TF1.cxx:1559
Ssiz_t Length() const
Definition: TString.h:405
Double_t * fParMax
Definition: TF1Data.h:49
int Status() const
return the Error Status of the last Integral calculation
User Class for performing numerical integration of a function in one dimension.
Definition: Integrator.h:94
Double_t fMinimum
Definition: TF1.h:249
static void RejectPoint(Bool_t reject=kTRUE)
Static function to set the global flag to reject points the fgRejectPoint global flag is tested by al...
Definition: TF1.cxx:3667
void GetParameters(TFitEditor::FuncParams_t &pars, TF1 *func)
Stores the parameters of the given function into pars.
Definition: TFitEditor.cxx:263
Class adding two functions: c1*f1+c2*f2.
Definition: TF1NormSum.h:19
static void InitStandardFunctions()
Create the basic function objects.
Definition: TF1.cxx:2486
virtual void SetMarkerStyle(Style_t mstyle=1)
Set the marker style.
Definition: TAttMarker.h:40
TAxis * GetYaxis()
Definition: TH1.h:317
float xmax
Definition: THbookFile.cxx:93
virtual void Update()
Called by functions such as SetRange, SetNpx, SetParameters to force the deletion of the associated h...
Definition: TF1.cxx:3632
virtual ~TF1()
TF1 default destructor.
Definition: TF1.cxx:938
virtual Int_t GetNpar() const
Definition: TFormula.h:238
virtual Double_t GetXmin() const
Definition: TF1.h:546
R__EXTERN TRandom * gRandom
Definition: TRandom.h:62
1-D histogram with a double per channel (see TH1 documentation)}
Definition: TH1.h:614
TString fFormula
Definition: TFormula.h:138
TString fName
Definition: TNamed.h:32
virtual TH1 * DoCreateHistogram(Double_t xmin, Double_t xmax, Bool_t recreate=kFALSE)
Create histogram with bin content equal to function value computed at the bin center This histogram w...
Definition: TF1.cxx:3039
Bool_t IsValid() const
Return true if the method call has been properly initialized and is usable.
virtual void SetMarkerSize(Size_t msize=1)
Set the marker size.
Definition: TAttMarker.h:41
static IntegrationMultiDim::Type DefaultIntegratorType()
TGraphErrors * gr
Definition: legend1.C:25
REAL epsilon
Definition: triangle.c:617
constexpr Double_t E()
Base of natural log: .
Definition: TMath.h:97
#define h(i)
Definition: RSha256.hxx:106
Int_t GetNpar() const
virtual Double_t GetMinimum(Double_t minval=-FLT_MAX) const
Return minimum value larger than minval of bins in the range, unless the value has been overridden by...
Definition: TH1.cxx:8069
void DefineNSUMTerm(TObjArray *newFuncs, TObjArray *coeffNames, TString &fullFormula, TString &formula, int termStart, int termEnd, Double_t xmin, Double_t xmax)
Helper functions for NSUM parsing.
Definition: TF1.cxx:867
void InitWithPrototype(TClass *cl, const char *method, const char *proto, Bool_t objectIsConst=kFALSE, ROOT::EFunctionMatchMode mode=ROOT::kConversionMatch)
Initialize the method invocation environment.
static void CalcGaussLegendreSamplingPoints(Int_t num, Double_t *x, Double_t *w, Double_t eps=3.0e-11)
Type safe interface (static method) The number of sampling points are taken from the TGraph...
Definition: TF1.cxx:3810
const Bool_t kFALSE
Definition: RtypesCore.h:88
bool R__SetClonesArrayTF1Updater(TF1Updater_t func)
PyObject * fType
Class for finding the root of a one dimensional function using the Brent algorithm.
virtual Color_t GetLineColor() const
Return the line color.
Definition: TAttLine.h:33
virtual void SetMaximum(Double_t maximum=-1111)
Set the maximum value along Y for this function In case the function is already drawn, set also the maximum in the helper histogram.
Definition: TF1.cxx:3410
virtual Int_t FindBin(Double_t x)
Find bin number corresponding to abscissa x.
Definition: TAxis.cxx:290
bool GetVectorizedOption(Option_t *opt)
Definition: TF1.cxx:680
virtual void SetTitle(const char *title="")
Set function title if title has the form "fffffff;xxxx;yyyy", it is assumed that the function title i...
Definition: TF1.cxx:3565
static unsigned int total
long Long_t
Definition: RtypesCore.h:50
static const std::string & DefaultMinimizerAlgo()
int R__RegisterTF1UpdaterTrigger
Definition: TF1.cxx:144
Double_t GetChisquare() const
Definition: TF1.h:438
virtual Double_t Eval(Double_t x, Double_t y=0, Double_t z=0, Double_t t=0) const
Evaluate this function.
Definition: TF1.cxx:1446
void SetMaxFunctionCalls(unsigned int maxfcn)
set maximum of function calls
Definition: Minimizer.h:448
class containg the result of the fit and all the related information (fitted parameter values...
Definition: FitResult.h:47
double Error() const
return integration error
Int_t fNDF
Definition: TF1.h:247
static const double x1[5]
#define ClassImp(name)
Definition: Rtypes.h:365
void SetFunction(const IGenFunction &)
Set integration function (flag control if function must be copied inside).
double Double_t
Definition: RtypesCore.h:55
std::vector< Double_t > fSave
Definition: TF1.h:254
void SetTolerance(double tol)
set the tolerance
Definition: Minimizer.h:454
std::vector< Double_t > fParErrors
Definition: TF1.h:251
Double_t fXmin
Definition: TF1.h:240
void SetParName(Int_t iparam, const char *name)
Definition: TF1.h:117
TF1FunctorPointer * fFunctor
Definition: TF1.h:264
Int_t GetParNumber(const char *name) const
Returns the parameter number given a name not very efficient but list of parameters is typically smal...
Definition: TF1.cxx:3830
Double_t AnalyticalIntegral(TF1 *f, Double_t a, Double_t b)
Double_t y[n]
Definition: legend1.C:17
std::vector< Double_t > fGamma
Array beta. is approximated by x = alpha +beta*r *gamma*r**2.
Definition: TF1.h:258
virtual Double_t GetX(Double_t y, Double_t xmin=0, Double_t xmax=0, Double_t epsilon=1.E-10, Int_t maxiter=100, Bool_t logx=false) const
Returns the X value corresponding to the function value fy for (xmin<x<xmax).
Definition: TF1.cxx:1860
virtual Color_t GetFillColor() const
Return the fill area color.
Definition: TAttFill.h:30
virtual Double_t GetXmax() const
Definition: TF1.h:550
void Streamer(TBuffer &b, Int_t version, UInt_t start, UInt_t count, const TClass *onfile_class=0)
Stream a class object.
Definition: TF1Data_v5.cxx:59
Bool_t Contains(const char *pat, ECaseCompare cmp=kExact) const
Definition: TString.h:619
The TH1 histogram class.
Definition: TH1.h:56
you should not use this method at all Int_t Int_t Double_t Double_t Double_t e
Definition: TRolke.cxx:630
static std::string GetName(IntegrationOneDim::Type)
static function to get a string from the enumeration
Definition: Integrator.cxx:66
#define R__LOCKGUARD(mutex)
virtual double DoEval(double x) const =0
implementation of the evaluation function. Must be implemented by derived classes ...
virtual void SavePrimitive(std::ostream &out, Option_t *option="")
Save primitive as a C++ statement(s) on output stream out.
Definition: TF1.cxx:3203
virtual Double_t Uniform(Double_t x1=1)
Returns a uniform deviate on the interval (0, x1).
Definition: TRandom.cxx:635
std::vector< Double_t > fAlpha
Integral of function binned on fNpx bins.
Definition: TF1.h:256
EAddToList
Definition: TF1.h:217
Double_t fXmin
Definition: TF1Data.h:39
virtual void SetLineStyle(Style_t lstyle)
Set the line style.
Definition: TAttLine.h:42
Template class to wrap any C++ callable object implementing operator() (const double * x) in a multi-...
Array of doubles (64 bits per element).
Definition: TArrayD.h:27
Namespace for new Math classes and functions.
Bool_t IsNull() const
Definition: TString.h:402
virtual void InitArgs(const Double_t *x, const Double_t *params)
Initialize parameters addresses.
Definition: TF1.cxx:2471
std::vector< Double_t > fBeta
Array alpha. for each bin in x the deconvolution r of fIntegral.
Definition: TF1.h:257
virtual TObject * Clone(const char *newname="") const
Make a clone of an object using the Streamer facility.
Definition: TNamed.cxx:74
void SetParNames(const char *name0="p0", const char *name1="p1", const char *name2="p2", const char *name3="p3", const char *name4="p4", const char *name5="p5", const char *name6="p6", const char *name7="p7", const char *name8="p8", const char *name9="p9", const char *name10="p10")
Set parameter names.
Definition: TF1.cxx:3862
virtual Double_t GetSave(const Double_t *x)
Get value corresponding to X in array of fSave values.
Definition: TF1.cxx:2343
bool IsEmpty() const
True if a fit result does not exist (even invalid) with parameter values.
Definition: FitResult.h:117
Mother of all ROOT objects.
Definition: TObject.h:37
you should not use this method at all Int_t Int_t z
Definition: TRolke.cxx:630
Double_t * fParMin
Definition: TF1Data.h:48
X-axis in log scale.
Definition: TH1.h:163
TMethodCall * fMethodCall
Pointer to histogram used for visualisation.
Definition: TF1.h:261
typedef void((*Func_t)())
Double_t fChisquare
Definition: TF1Data.h:46
virtual Int_t GetNpar() const
Definition: TF1.h:475
virtual Double_t IntegralOneDim(Double_t a, Double_t b, Double_t epsrel, Double_t epsabs, Double_t &err)
Return Integral of function between a and b using the given parameter values and relative and absolut...
Definition: TF1.cxx:2609
virtual void GetParLimits(Int_t ipar, Double_t &parmin, Double_t &parmax) const
Return limits for parameter ipar.
Definition: TF1.cxx:1937
virtual bool SetVariable(unsigned int ivar, const std::string &name, double val, double step)=0
set a new free variable
virtual void Copy(TObject &named) const
Copy this to obj.
Definition: TNamed.cxx:94
User class for performing multidimensional integration.
static Bool_t RejectedPoint()
See TF1::RejectPoint above.
Definition: TF1.cxx:3676
Style_t GetFuncStyle() const
Definition: TStyle.h:209
auto * l
Definition: textangle.C:4
Double_t fMaximum
Definition: TF1Data.h:51
virtual Double_t GetMinMaxNDim(Double_t *x, Bool_t findmax, Double_t epsilon=0, Int_t maxiter=0) const
Find the minimum of a function of whatever dimension.
Definition: TF1.cxx:1720
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:3500
void Execute(const char *, const char *, int *=0)
Execute method on this object with the given parameter string, e.g.
Definition: TMethodCall.h:64
virtual Double_t GetParameter(Int_t ipar) const
Definition: TF1.h:506
virtual void GetRange(Double_t *xmin, Double_t *xmax) const
Return range of a generic N-D function.
Definition: TF1.cxx:2280
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:831
Short_t Max(Short_t a, Short_t b)
Definition: TMathBase.h:212
1-Dim function class
Definition: TF1.h:210
void MakeZombie()
Definition: TObject.h:49
virtual void SetSavedPoint(Int_t point, Double_t value)
Restore value of function saved at point.
Definition: TF1.cxx:3549
A TGraph is an object made of two arrays X and Y with npoints each.
Definition: TGraph.h:41
TObject * Clone(const char *newname=0) const
Make a complete copy of the underlying object.
Definition: TH1.cxx:2663
virtual void SetFitResult(const ROOT::Fit::FitResult &result, const Int_t *indpar=0)
Set the result from the fit parameter values, errors, chi2, etc...
Definition: TF1.cxx:3371
TF1 * f1
Definition: legend1.C:11
you should not use this method at all Int_t Int_t Double_t Double_t Double_t Int_t Double_t Double_t Double_t Double_t b
Definition: TRolke.cxx:630
virtual Double_t * GetParameters() const
Definition: TFormula.h:243
const Double_t * GetParameters() const
Definition: TF1.h:82
TF1Parameters * fParams
Definition: TF1.h:266
static std::atomic< Bool_t > fgAddToGlobList
Definition: TF1.h:302
virtual void Print(Option_t *option="") const
Print TNamed name and title.
Definition: TF1.cxx:2884
#define snprintf
Definition: civetweb.c:1540
Color_t GetFuncColor() const
Definition: TStyle.h:208
#define gPad
Definition: TVirtualPad.h:287
static void AbsValue(Bool_t reject=kTRUE)
Static function: set the fgAbsValue flag.
Definition: TF1.cxx:973
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:3435
R__EXTERN Int_t gDebug
Definition: Rtypes.h:91
virtual IBaseFunctionOneDim * Clone() const =0
Clone a function.
virtual Double_t * GetParameters() const
Definition: TF1.h:514
void Add(TObject *obj)
Definition: TObjArray.h:74
virtual Double_t IntegralError(Double_t a, Double_t b, const Double_t *params=0, const Double_t *covmat=0, Double_t epsilon=1.E-2)
Return Error on Integral of a parametric function between a and b due to the parameter uncertainties...
Definition: TF1.cxx:2706
virtual TH1 * GetHistogram() const
Return a pointer to the histogram used to visualise the function.
Definition: TF1.cxx:1581
TF1 & operator=(const TF1 &rhs)
Operator =.
Definition: TF1.cxx:926
virtual void SetParameter(Int_t param, Double_t value)
Definition: TF1.h:628
virtual void SetTitle(const char *title)
See GetStatOverflows for more information.
Definition: TH1.cxx:6323
bool SetFunction(const ROOT::Math::IGenFunction &f, double xlow, double xup)
Sets the function for the rest of the algorithms.
void SetNpx(int npx)
Set the number of point used to bracket root using a grid.
virtual Color_t GetMarkerColor() const
Return the marker color.
Definition: TAttMarker.h:31
ParamFunctorTempl< double > ParamFunctor
Definition: ParamFunctor.h:388
Double_t * fSave
Definition: TF1Data.h:50
Functor1D class for one-dimensional functions.
Definition: Functor.h:493
virtual Style_t GetFillStyle() const
Return the fill area style.
Definition: TAttFill.h:31
Double_t Sqrt(Double_t x)
Definition: TMath.h:681
virtual TString GetExpFormula(Option_t *option="") const
Reconstruct the formula expression from the internal TFormula member variables.
std::vector< Double_t > fParMax
Definition: TF1.h:253
virtual Double_t CentralMoment(Double_t n, Double_t a, Double_t b, const Double_t *params=0, Double_t epsilon=0.000001)
Return nth central moment of function between a and b (i.e the n-th moment around the mean value) ...
Definition: TF1.cxx:3723
float * q
Definition: THbookFile.cxx:87
void Print(Option_t *option="") const
Print TNamed name and title.
virtual Double_t EvalPar(const Double_t *x, const Double_t *params=0)
Evaluate function with given coordinates and parameters.
Definition: TF1.cxx:1475
virtual void SetTitle(const char *title="")
Set the title of the TNamed.
Definition: TNamed.cxx:164
Bool_t IsValid() const
Definition: TFormula.h:233
virtual void Browse(TBrowser *b)
Browse.
Definition: TF1.cxx:982
double RelError() const
return relative error
const Bool_t kTRUE
Definition: RtypesCore.h:87
Double_t GetXmax() const
Definition: TAxis.h:134
virtual Int_t GetNumber() const
Definition: TF1.h:492
double Error() const
Returns the estimate of the absolute Error of the last derivative calculation.
std::vector< Double_t > fParMin
Definition: TF1.h:252
const Int_t n
Definition: legend1.C:16
Line Attributes class.
Definition: TAttLine.h:18
Double_t fXmax
Definition: TF1.h:241
Long64_t BinarySearch(Long64_t n, const T *array, T value)
Definition: TMathBase.h:278
static Double_t gErrorTF1
Definition: TF1.cxx:56
double Chi2() const
Chi2 fit value in case of likelihood must be computed ?
Definition: FitResult.h:160
User class for calculating the derivatives of a function.
char name[80]
Definition: TGX11.cxx:109
Class for adaptive quadrature integration in multi-dimensions using rectangular regions.
virtual void Warning(const char *method, const char *msgfmt,...) const
Issue warning message.
Definition: TObject.cxx:873
TAxis * GetXaxis()
Get the behaviour adopted by the object about the statoverflows. See EStatOverflows for more informat...
Definition: TH1.h:316
virtual Int_t GetParNumber(const char *name) const
Definition: TF1.h:527
virtual const char * GetTitle() const
Returns title of object.
Definition: TNamed.h:48
Double_t fMinimum
Definition: TF1Data.h:52
void GetWeightVectors(double *x, double *w) const
Returns the arrays x and w containing the abscissa and weight of the Gauss-Legendre n-point quadratur...
TAxis * GetZaxis() const
Get z axis of the function. (In case this object is a TF2 or TF3)
Definition: TF1.cxx:2421
void Copy(TAttFill &attfill) const
Copy this fill attributes to a new TAttFill.
Definition: TAttFill.cxx:202
virtual void Paint(Option_t *option="")
Paint this function with its current attributes.
Definition: TF1.cxx:2940
const char * Data() const
Definition: TString.h:364
virtual TObject * DrawDerivative(Option_t *option="al")
Draw derivative of this function.
Definition: TF1.cxx:1386