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