Logo ROOT   6.30.04
Reference Guide
 All Namespaces Files Pages
TF1.h
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 // ---------------------------------- F1.h
12 
13 #ifndef ROOT_TF1
14 #define ROOT_TF1
15 
16 //////////////////////////////////////////////////////////////////////////
17 // //
18 // TF1 //
19 // //
20 // The Parametric 1-D function //
21 // //
22 //////////////////////////////////////////////////////////////////////////
23 
24 #include "RConfigure.h"
25 #include <functional>
26 #include <cassert>
27 #include "TFormula.h"
28 #include "TAttLine.h"
29 #include "TAttFill.h"
30 #include "TAttMarker.h"
31 #include "TROOT.h"
32 #include "TF1AbsComposition.h"
33 #include "TMath.h"
34 #include "Math/Types.h"
35 #include "Math/ParamFunctor.h"
36 
37 class TF1;
38 class TH1;
39 class TAxis;
40 class TMethodCall;
41 
42 namespace ROOT {
43  namespace Fit {
44  class FitResult;
45  }
46 }
47 
48 class TF1Parameters {
49 public:
50  TF1Parameters() {} // needed for the I/O
51  TF1Parameters(Int_t npar) :
52  fParameters(std::vector<Double_t>(npar)),
53  fParNames(std::vector<std::string>(npar))
54  {
55  for (int i = 0; i < npar; ++i) {
56  fParNames[i] = std::string(TString::Format("p%d", i).Data());
57  }
58  }
59  // copy constructor
60  TF1Parameters(const TF1Parameters &rhs) :
61  fParameters(rhs.fParameters),
62  fParNames(rhs.fParNames)
63  {}
64  // assignment
65  TF1Parameters &operator=(const TF1Parameters &rhs)
66  {
67  if (&rhs == this) return *this;
68  fParameters = rhs.fParameters;
69  fParNames = rhs.fParNames;
70  return *this;
71  }
72  virtual ~TF1Parameters() {}
73 
74  // getter methods
75  Double_t GetParameter(Int_t iparam) const
76  {
77  return (CheckIndex(iparam)) ? fParameters[iparam] : 0;
78  }
79  Double_t GetParameter(const char *name) const
80  {
81  return GetParameter(GetParNumber(name));
82  }
83  const Double_t *GetParameters() const
84  {
85  return fParameters.data();
86  }
87  const std::vector<double> &ParamsVec() const
88  {
89  return fParameters;
90  }
91 
92  Int_t GetParNumber(const char *name) const;
93 
94  const char *GetParName(Int_t iparam) const
95  {
96  return (CheckIndex(iparam)) ? fParNames[iparam].c_str() : "";
97  }
98 
99 
100  // setter methods
101  void SetParameter(Int_t iparam, Double_t value)
102  {
103  if (!CheckIndex(iparam)) return;
104  fParameters[iparam] = value;
105  }
106  void SetParameters(const Double_t *params)
107  {
108  std::copy(params, params + fParameters.size(), fParameters.begin());
109  }
110  void SetParameters(Double_t p0, Double_t p1, Double_t p2 = 0, Double_t p3 = 0, Double_t p4 = 0,
111  Double_t p5 = 0, Double_t p6 = 0, Double_t p7 = 0, Double_t p8 = 0,
112  Double_t p9 = 0, Double_t p10 = 0);
113 
114  void SetParameter(const char *name, Double_t value)
115  {
116  SetParameter(GetParNumber(name), value);
117  }
118  void SetParName(Int_t iparam, const char *name)
119  {
120  if (!CheckIndex(iparam)) return;
121  fParNames[iparam] = std::string(name);
122  }
123  void SetParNames(const char *name0 = "p0", const char *name1 = "p1", const char *name2 = "p2",
124  const char *name3 = "p3", const char *name4 = "p4", const char *name5 = "p5",
125  const char *name6 = "p6", const char *name7 = "p7", const char *name8 = "p8",
126  const char *name9 = "p9", const char *name10 = "p10");
127 
128 
129 
130  ClassDef(TF1Parameters, 1) // The Parameters of a parameteric function
131 private:
132 
133  bool CheckIndex(Int_t i) const
134  {
135  return (i >= 0 && i < int(fParameters.size()));
136  }
137 
138  std::vector<Double_t> fParameters; // parameter values
139  std::vector<std::string> fParNames; // parameter names
140 };
141 
142 namespace ROOT {
143  namespace Internal {
144  /// Internal class used by TF1 for defining
145  /// template specialization for different TF1 constructors
146  template<class Func>
147  struct TF1Builder {
148  static void Build(TF1 *f, Func func);
149  };
150 
151  template<class Func>
152  struct TF1Builder<Func *> {
153  static void Build(TF1 *f, Func *func);
154  };
155 
156  // Internal class used by TF1 for obtaining the type from a functor
157  // out of the set of valid operator() signatures.
158  template<typename T>
159  struct GetFunctorType {
160  };
161 
162  template<typename F, typename T>
163  struct GetFunctorType<T(F::*)(const T *, const double *)> {
164  using type = T;
165  };
166 
167  template<typename F, typename T>
168  struct GetFunctorType<T(F::*)(const T *, const double *) const> {
169  using type = T;
170  };
171 
172  template<typename F, typename T>
173  struct GetFunctorType<T(F::*)(T *, double *)> {
174  using type = T;
175  };
176 
177  template<typename F, typename T>
178  struct GetFunctorType<T(F::*)(T *, double *) const> {
179  using type = T;
180  };
181 
182  // Internal class used by TF1 to get the right operator() signature
183  // from a Functor with several ones.
184  template<typename T, typename F>
185  auto GetTheRightOp(T(F::*opPtr)(const T *, const double *)) -> decltype(opPtr)
186  {
187  return opPtr;
188  }
189 
190  template<typename T, typename F>
191  auto GetTheRightOp(T(F::*opPtr)(const T *, const double *) const) -> decltype(opPtr)
192  {
193  return opPtr;
194  }
195 
196  template<typename T, typename F>
197  auto GetTheRightOp(T(F::*opPtr)(T *, double *)) -> decltype(opPtr)
198  {
199  return opPtr;
200  }
201 
202  template<typename T, typename F>
203  auto GetTheRightOp(T(F::*opPtr)(T *, double *) const) -> decltype(opPtr)
204  {
205  return opPtr;
206  }
207  }
208 }
209 
210 
211 class TF1 : public TNamed, public TAttLine, public TAttFill, public TAttMarker {
212 
213  template<class Func>
214  friend struct ROOT::Internal::TF1Builder;
215 
216 public:
217  // Add to list behavior
218  enum class EAddToList {
219  kDefault,
220  kAdd,
221  kNo
222  };
223 
224 protected:
225 
226  struct TF1FunctorPointer {
227  virtual ~TF1FunctorPointer() {}
228  virtual TF1FunctorPointer * Clone() const = 0;
229  };
230 
231 
232  enum EFType {
233  kFormula = 0, // formula functions which can be stored,
234  kPtrScalarFreeFcn, // pointer to scalar free function,
235  kInterpreted, // interpreted functions constructed by name,
236  kTemplVec, // vectorized free functions or TemplScalar functors evaluating on vectorized parameters,
237  kTemplScalar, // TemplScalar functors evaluating on scalar parameters
238  kCompositionFcn
239  }; // formula based on composition class (e.g. NSUM, CONV)
240 
241  Double_t fXmin{-1111}; //Lower bounds for the range
242  Double_t fXmax{-1111}; //Upper bounds for the range
243  Int_t fNpar{}; //Number of parameters
244  Int_t fNdim{}; //Function dimension
245  Int_t fNpx{100}; //Number of points used for the graphical representation
246  EFType fType{EFType::kTemplScalar};
247  Int_t fNpfits{}; //Number of points used in the fit
248  Int_t fNDF{}; //Number of degrees of freedom in the fit
249  Double_t fChisquare{}; //Function fit chisquare
250  Double_t fMinimum{-1111}; //Minimum value for plotting
251  Double_t fMaximum{-1111}; //Maximum value for plotting
252  std::vector<Double_t> fParErrors; //Array of errors of the fNpar parameters
253  std::vector<Double_t> fParMin; //Array of lower limits of the fNpar parameters
254  std::vector<Double_t> fParMax; //Array of upper limits of the fNpar parameters
255  std::vector<Double_t> fSave; //Array of fNsave function values
256  std::vector<Double_t> fIntegral; //!Integral of function binned on fNpx bins
257  std::vector<Double_t> fAlpha; //!Array alpha. for each bin in x the deconvolution r of fIntegral
258  std::vector<Double_t> fBeta; //!Array beta. is approximated by x = alpha +beta*r *gamma*r**2
259  std::vector<Double_t> fGamma; //!Array gamma.
260  TObject *fParent{nullptr}; //!Parent object hooking this function (if one)
261  TH1 *fHistogram{nullptr}; //!Pointer to histogram used for visualisation
262  TMethodCall *fMethodCall{nullptr}; //!Pointer to MethodCall in case of interpreted function
263  Bool_t fNormalized{false}; //Normalization option (false by default)
264  Double_t fNormIntegral{}; //Integral of the function before being normalized
265  TF1FunctorPointer *fFunctor{nullptr}; //! Functor object to wrap any C++ callable object
266  TFormula *fFormula{nullptr}; //Pointer to TFormula in case when user define formula
267  TF1Parameters *fParams{nullptr}; //Pointer to Function parameters object (exists only for not-formula functions)
268  std::unique_ptr<TF1AbsComposition> fComposition; //! Pointer to composition (NSUM or CONV)
269  TF1AbsComposition *fComposition_ptr{nullptr}; // saved pointer (unique_ptr is transient)
270 
271  /// General constructor for TF1. Most of the other constructors delegate on it
272  TF1(EFType functionType, const char *name, Double_t xmin, Double_t xmax, Int_t npar, Int_t ndim, EAddToList addToGlobList, TF1Parameters *params = nullptr, TF1FunctorPointer * functor = nullptr):
273  TNamed(name, name), TAttLine(), TAttFill(), TAttMarker(), fXmin(xmin), fXmax(xmax), fNpar(npar), fNdim(ndim),
274  fType(functionType), fParErrors(npar), fParMin(npar), fParMax(npar), fFunctor(functor), fParams(params)
275  {
276  DoInitialize(addToGlobList);
277  };
278 
279 private:
280  // NSUM parsing helper functions
281  void DefineNSUMTerm(TObjArray *newFuncs, TObjArray *coeffNames,
282  TString &fullFormula,
283  TString &formula, int termStart, int termEnd,
284  Double_t xmin, Double_t xmax);
285  int TermCoeffLength(TString &term);
286 
287 protected:
288 
289  template <class T>
290  struct TF1FunctorPointerImpl: TF1FunctorPointer {
291  TF1FunctorPointerImpl(const ROOT::Math::ParamFunctorTempl<T> &func): fImpl(func) {};
292  TF1FunctorPointerImpl(const std::function<T(const T *f, const Double_t *param)> &func) : fImpl(func){};
293  virtual ~TF1FunctorPointerImpl() {}
294  virtual TF1FunctorPointer * Clone() const { return new TF1FunctorPointerImpl<T>(fImpl); }
295  ROOT::Math::ParamFunctorTempl<T> fImpl;
296  };
297 
298 
299 
300 
301  static std::atomic<Bool_t> fgAbsValue; //use absolute value of function when computing integral
302  static Bool_t fgRejectPoint; //True if point must be rejected in a fit
303  static std::atomic<Bool_t> fgAddToGlobList; //True if we want to register the function in the global list
304  static TF1 *fgCurrent; //pointer to current function being processed
305 
306 
307  //void CreateFromFunctor(const char *name, Int_t npar, Int_t ndim = 1);
308  void DoInitialize(EAddToList addToGlobList);
309 
310  void IntegrateForNormalization();
311 
312  virtual Double_t GetMinMaxNDim(Double_t *x , Bool_t findmax, Double_t epsilon = 0, Int_t maxiter = 0) const;
313  virtual void GetRange(Double_t *xmin, Double_t *xmax) const;
314  virtual TH1 *DoCreateHistogram(Double_t xmin, Double_t xmax, Bool_t recreate = kFALSE);
315 
316 public:
317 
318  // TF1 status bits
319  enum EStatusBits {
320  kNotGlobal = BIT(10), // don't register in global list of functions
321  kNotDraw = BIT(9) // don't draw the function when in a TH1
322  };
323 
324  TF1();
325  TF1(const char *name, const char *formula, Double_t xmin = 0, Double_t xmax = 1, EAddToList addToGlobList = EAddToList::kDefault, bool vectorize = false);
326  TF1(const char *name, const char *formula, Double_t xmin, Double_t xmax, Option_t * option); // same as above but using a string for option
327  TF1(const char *name, Double_t xmin, Double_t xmax, Int_t npar, Int_t ndim = 1, EAddToList addToGlobList = EAddToList::kDefault);
328  TF1(const char *name, Double_t (*fcn)(Double_t *, Double_t *), Double_t xmin = 0, Double_t xmax = 1, Int_t npar = 0, Int_t ndim = 1, EAddToList addToGlobList = EAddToList::kDefault);
329  TF1(const char *name, Double_t (*fcn)(const Double_t *, const Double_t *), Double_t xmin = 0, Double_t xmax = 1, Int_t npar = 0, Int_t ndim = 1, EAddToList addToGlobList = EAddToList::kDefault);
330 
331  template <class T>
332  TF1(const char *name, std::function<T(const T *data, const Double_t *param)> &fcn, Double_t xmin = 0, Double_t xmax = 1, Int_t npar = 0, Int_t ndim = 1, EAddToList addToGlobList = EAddToList::kDefault):
333  TF1(EFType::kTemplScalar, name, xmin, xmax, npar, ndim, addToGlobList, new TF1Parameters(npar), new TF1FunctorPointerImpl<T>(fcn))
334  {
335  fType = std::is_same<T, double>::value ? TF1::EFType::kTemplScalar : TF1::EFType::kTemplVec;
336  }
337 
338  ////////////////////////////////////////////////////////////////////////////////
339  /// Constructor using a pointer to function.
340  ///
341  /// \param npar is the number of free parameters used by the function
342  ///
343  /// This constructor creates a function of type C when invoked
344  /// with the normal C++ compiler.
345  ///
346  ///
347  /// WARNING! A function created with this constructor cannot be Cloned
348 
349 
350  template <class T>
351  TF1(const char *name, T(*fcn)(const T *, const Double_t *), Double_t xmin = 0, Double_t xmax = 1, Int_t npar = 0, Int_t ndim = 1, EAddToList addToGlobList = EAddToList::kDefault):
352  TF1(EFType::kTemplVec, name, xmin, xmax, npar, ndim, addToGlobList, new TF1Parameters(npar), new TF1FunctorPointerImpl<T>(fcn))
353  {}
354 
355  // Constructors using functors (compiled mode only)
356  TF1(const char *name, ROOT::Math::ParamFunctor f, Double_t xmin = 0, Double_t xmax = 1, Int_t npar = 0, Int_t ndim = 1, EAddToList addToGlobList = EAddToList::kDefault);
357 
358  // Template constructors from any C++ callable object, defining the operator() (double * , double *)
359  // and returning a double.
360  // The class name is not needed when using compile code, while it is required when using
361  // interpreted code via the specialized constructor with void *.
362  // An instance of the C++ function class or its pointer can both be used. The former is reccomended when using
363  // C++ compiled code, but if CINT compatibility is needed, then a pointer to the function class must be used.
364  // xmin and xmax specify the plotting range, npar is the number of parameters.
365  // See the tutorial math/exampleFunctor.C for an example of using this constructor
366  template <typename Func>
367  TF1(const char *name, Func f, Double_t xmin, Double_t xmax, Int_t npar, Int_t ndim = 1, EAddToList addToGlobList = EAddToList::kDefault) :
368  TF1(EFType::kTemplScalar, name, xmin, xmax, npar, ndim, addToGlobList)
369  {
370  //actual fType set in TF1Builder
371  ROOT::Internal::TF1Builder<Func>::Build(this, f);
372  }
373 
374  // backward compatible interface
375  template <typename Func>
376  TF1(const char *name, Func f, Double_t xmin, Double_t xmax, Int_t npar, const char *, EAddToList addToGlobList = EAddToList::kDefault) :
377  TF1(EFType::kTemplScalar, name, xmin, xmax, npar, 1, addToGlobList, new TF1Parameters(npar))
378  {
379  ROOT::Internal::TF1Builder<Func>::Build(this, f);
380  }
381 
382 
383  // Template constructors from a pointer to any C++ class of type PtrObj with a specific member function of type
384  // MemFn.
385  // The member function must have the signature of (double * , double *) and returning a double.
386  // The class name and the method name are not needed when using compile code
387  // (the member function pointer is used in this case), while they are required when using interpreted
388  // code via the specialized constructor with void *.
389  // xmin and xmax specify the plotting range, npar is the number of parameters.
390  // See the tutorial math/exampleFunctor.C for an example of using this constructor
391  template <class PtrObj, typename MemFn>
392  TF1(const char *name, const PtrObj &p, MemFn memFn, Double_t xmin, Double_t xmax, Int_t npar, Int_t ndim = 1, EAddToList addToGlobList = EAddToList::kDefault) :
393  TF1(EFType::kTemplScalar, name, xmin, xmax, npar, ndim, addToGlobList, new TF1Parameters(npar), new TF1FunctorPointerImpl<double>(ROOT::Math::ParamFunctor(p, memFn)))
394  {}
395 
396  // backward compatible interface
397  template <class PtrObj, typename MemFn>
398  TF1(const char *name, const PtrObj &p, MemFn memFn, Double_t xmin, Double_t xmax, Int_t npar, const char *, const char *, EAddToList addToGlobList = EAddToList::kDefault) :
399  TF1(EFType::kTemplScalar, name, xmin, xmax, npar, 1, addToGlobList, new TF1Parameters(npar), new TF1FunctorPointerImpl<double>(ROOT::Math::ParamFunctor(p, memFn)))
400  {}
401 
402  TF1(const TF1 &f1);
403  TF1 &operator=(const TF1 &rhs);
404  virtual ~TF1();
405  virtual void AddParameter(const TString &name, Double_t value)
406  {
407  if (fFormula) fFormula->AddParameter(name, value);
408  }
409  // virtual void AddParameters(const pair<TString,Double_t> *pairs, Int_t size) { fFormula->AddParameters(pairs,size); }
410  // virtual void AddVariable(const TString &name, Double_t value = 0) { if (fFormula) fFormula->AddVariable(name,value); }
411  // virtual void AddVariables(const TString *vars, Int_t size) { if (fFormula) fFormula->AddVariables(vars,size); }
412  virtual Bool_t AddToGlobalList(Bool_t on = kTRUE);
413  static Bool_t DefaultAddToGlobalList(Bool_t on = kTRUE);
414  virtual void Browse(TBrowser *b);
415  virtual void Copy(TObject &f1) const;
416  virtual Double_t Derivative(Double_t x, Double_t *params = 0, Double_t epsilon = 0.001) const;
417  virtual Double_t Derivative2(Double_t x, Double_t *params = 0, Double_t epsilon = 0.001) const;
418  virtual Double_t Derivative3(Double_t x, Double_t *params = 0, Double_t epsilon = 0.001) const;
419  static Double_t DerivativeError();
420  virtual Int_t DistancetoPrimitive(Int_t px, Int_t py);
421  virtual void Draw(Option_t *option = "");
422  virtual TF1 *DrawCopy(Option_t *option = "") const;
423  virtual TObject *DrawDerivative(Option_t *option = "al"); // *MENU*
424  virtual TObject *DrawIntegral(Option_t *option = "al"); // *MENU*
425  virtual void DrawF1(Double_t xmin, Double_t xmax, Option_t *option = "");
426  virtual Double_t Eval(Double_t x, Double_t y = 0, Double_t z = 0, Double_t t = 0) const;
427  //template <class T> T Eval(T x, T y = 0, T z = 0, T t = 0) const;
428  virtual Double_t EvalPar(const Double_t *x, const Double_t *params = 0);
429  template <class T> T EvalPar(const T *x, const Double_t *params = 0);
430  virtual Double_t operator()(Double_t x, Double_t y = 0, Double_t z = 0, Double_t t = 0) const;
431  template <class T> T operator()(const T *x, const Double_t *params = nullptr);
432  virtual void ExecuteEvent(Int_t event, Int_t px, Int_t py);
433  virtual void FixParameter(Int_t ipar, Double_t value);
434  bool IsVectorized()
435  {
436  return (fType == EFType::kTemplVec) || (fType == EFType::kFormula && fFormula && fFormula->IsVectorized());
437  }
438  Double_t GetChisquare() const
439  {
440  return fChisquare;
441  }
442  virtual TH1 *GetHistogram() const;
443  virtual TH1 *CreateHistogram()
444  {
445  return DoCreateHistogram(fXmin, fXmax);
446  }
447  virtual TFormula *GetFormula()
448  {
449  return fFormula;
450  }
451  virtual const TFormula *GetFormula() const
452  {
453  return fFormula;
454  }
455  virtual TString GetExpFormula(Option_t *option = "") const
456  {
457  return (fFormula) ? fFormula->GetExpFormula(option) : "";
458  }
459  virtual const TObject *GetLinearPart(Int_t i) const
460  {
461  return (fFormula) ? fFormula->GetLinearPart(i) : nullptr;
462  }
463  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;
464  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;
465  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;
466  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;
467  virtual Double_t GetMaximumStored() const
468  {
469  return fMaximum;
470  }
471  virtual Double_t GetMinimumStored() const
472  {
473  return fMinimum;
474  }
475  virtual Int_t GetNpar() const
476  {
477  return fNpar;
478  }
479  virtual Int_t GetNdim() const
480  {
481  return fNdim;
482  }
483  virtual Int_t GetNDF() const;
484  virtual Int_t GetNpx() const
485  {
486  return fNpx;
487  }
488  TMethodCall *GetMethodCall() const
489  {
490  return fMethodCall;
491  }
492  virtual Int_t GetNumber() const
493  {
494  return (fFormula) ? fFormula->GetNumber() : 0;
495  }
496  virtual Int_t GetNumberFreeParameters() const;
497  virtual Int_t GetNumberFitPoints() const
498  {
499  return fNpfits;
500  }
501  virtual char *GetObjectInfo(Int_t px, Int_t py) const;
502  TObject *GetParent() const
503  {
504  return fParent;
505  }
506  virtual Double_t GetParameter(Int_t ipar) const
507  {
508  return (fFormula) ? fFormula->GetParameter(ipar) : fParams->GetParameter(ipar);
509  }
510  virtual Double_t GetParameter(const TString &name) const
511  {
512  return (fFormula) ? fFormula->GetParameter(name) : fParams->GetParameter(name);
513  }
514  virtual Double_t *GetParameters() const
515  {
516  return (fFormula) ? fFormula->GetParameters() : const_cast<Double_t *>(fParams->GetParameters());
517  }
518  virtual void GetParameters(Double_t *params)
519  {
520  if (fFormula) fFormula->GetParameters(params);
521  else std::copy(fParams->ParamsVec().begin(), fParams->ParamsVec().end(), params);
522  }
523  virtual const char *GetParName(Int_t ipar) const
524  {
525  return (fFormula) ? fFormula->GetParName(ipar) : fParams->GetParName(ipar);
526  }
527  virtual Int_t GetParNumber(const char *name) const
528  {
529  return (fFormula) ? fFormula->GetParNumber(name) : fParams->GetParNumber(name);
530  }
531  virtual Double_t GetParError(Int_t ipar) const;
532  virtual const Double_t *GetParErrors() const
533  {
534  return fParErrors.data();
535  }
536  virtual void GetParLimits(Int_t ipar, Double_t &parmin, Double_t &parmax) const;
537  virtual Double_t GetProb() const;
538  virtual Int_t GetQuantiles(Int_t nprobSum, Double_t *q, const Double_t *probSum);
539  virtual Double_t GetRandom();
540  virtual Double_t GetRandom(Double_t xmin, Double_t xmax);
541  virtual void GetRange(Double_t &xmin, Double_t &xmax) const;
542  virtual void GetRange(Double_t &xmin, Double_t &ymin, Double_t &xmax, Double_t &ymax) const;
543  virtual void GetRange(Double_t &xmin, Double_t &ymin, Double_t &zmin, Double_t &xmax, Double_t &ymax, Double_t &zmax) const;
544  virtual Double_t GetSave(const Double_t *x);
545  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;
546  virtual Double_t GetXmin() const
547  {
548  return fXmin;
549  }
550  virtual Double_t GetXmax() const
551  {
552  return fXmax;
553  }
554  TAxis *GetXaxis() const ;
555  TAxis *GetYaxis() const ;
556  TAxis *GetZaxis() const ;
557  virtual Double_t GetVariable(const TString &name)
558  {
559  return (fFormula) ? fFormula->GetVariable(name) : 0;
560  }
561  virtual Double_t GradientPar(Int_t ipar, const Double_t *x, Double_t eps = 0.01);
562  template <class T>
563  T GradientPar(Int_t ipar, const T *x, Double_t eps = 0.01);
564  template <class T>
565  T GradientParTempl(Int_t ipar, const T *x, Double_t eps = 0.01);
566 
567  virtual void GradientPar(const Double_t *x, Double_t *grad, Double_t eps = 0.01);
568  template <class T>
569  void GradientPar(const T *x, T *grad, Double_t eps = 0.01);
570  template <class T>
571  void GradientParTempl(const T *x, T *grad, Double_t eps = 0.01);
572 
573  virtual void InitArgs(const Double_t *x, const Double_t *params);
574  static void InitStandardFunctions();
575  virtual Double_t Integral(Double_t a, Double_t b, Double_t epsrel = 1.e-12);
576  virtual Double_t IntegralOneDim(Double_t a, Double_t b, Double_t epsrel, Double_t epsabs, Double_t &err);
577  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);
578  virtual Double_t IntegralError(Int_t n, const Double_t *a, const Double_t *b, const Double_t *params = 0, const Double_t *covmat = 0, Double_t epsilon = 1.E-2);
579  // virtual Double_t IntegralFast(const TGraph *g, Double_t a, Double_t b, Double_t *params=0);
580  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);
581  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);
582  virtual Double_t IntegralMultiple(Int_t n, const Double_t *a, const Double_t *b, Int_t /*minpts*/, Int_t maxpts, Double_t epsrel, Double_t &relerr, Int_t &nfnevl, Int_t &ifail)
583  {
584  return IntegralMultiple(n, a, b, maxpts, epsrel, epsrel, relerr, nfnevl, ifail);
585  }
586  virtual Double_t IntegralMultiple(Int_t n, const Double_t *a, const Double_t *b, Double_t epsrel, Double_t &relerr);
587  virtual Bool_t IsEvalNormalized() const
588  {
589  return fNormalized;
590  }
591  /// return kTRUE if the point is inside the function range
592  virtual Bool_t IsInside(const Double_t *x) const
593  {
594  return !((x[0] < fXmin) || (x[0] > fXmax));
595  }
596  virtual Bool_t IsLinear() const
597  {
598  return (fFormula) ? fFormula->IsLinear() : false;
599  }
600  virtual Bool_t IsValid() const;
601  virtual void Print(Option_t *option = "") const;
602  virtual void Paint(Option_t *option = "");
603  virtual void ReleaseParameter(Int_t ipar);
604  virtual void Save(Double_t xmin, Double_t xmax, Double_t ymin, Double_t ymax, Double_t zmin, Double_t zmax);
605  virtual void SavePrimitive(std::ostream &out, Option_t *option = "");
606  virtual void SetChisquare(Double_t chi2)
607  {
608  fChisquare = chi2;
609  }
610  virtual void SetFitResult(const ROOT::Fit::FitResult &result, const Int_t *indpar = 0);
611  template <class PtrObj, typename MemFn>
612  void SetFunction(PtrObj &p, MemFn memFn);
613  template <typename Func>
614  void SetFunction(Func f);
615  virtual void SetMaximum(Double_t maximum = -1111); // *MENU*
616  virtual void SetMinimum(Double_t minimum = -1111); // *MENU*
617  virtual void SetNDF(Int_t ndf);
618  virtual void SetNumberFitPoints(Int_t npfits)
619  {
620  fNpfits = npfits;
621  }
622  virtual void SetNormalized(Bool_t flag)
623  {
624  fNormalized = flag;
625  Update();
626  }
627  virtual void SetNpx(Int_t npx = 100); // *MENU*
628  virtual void SetParameter(Int_t param, Double_t value)
629  {
630  (fFormula) ? fFormula->SetParameter(param, value) : fParams->SetParameter(param, value);
631  Update();
632  }
633  virtual void SetParameter(const TString &name, Double_t value)
634  {
635  (fFormula) ? fFormula->SetParameter(name, value) : fParams->SetParameter(name, value);
636  Update();
637  }
638  virtual void SetParameters(const Double_t *params)
639  {
640  (fFormula) ? fFormula->SetParameters(params) : fParams->SetParameters(params);
641  Update();
642  }
643  virtual void SetParameters(Double_t p0, Double_t p1, Double_t p2 = 0, Double_t p3 = 0, Double_t p4 = 0,
644  Double_t p5 = 0, Double_t p6 = 0, Double_t p7 = 0, Double_t p8 = 0,
645  Double_t p9 = 0, Double_t p10 = 0)
646  {
647  if (fFormula) fFormula->SetParameters(p0, p1, p2, p3, p4, p5, p6, p7, p8, p9, p10);
648  else fParams->SetParameters(p0, p1, p2, p3, p4, p5, p6, p7, p8, p9, p10);
649  Update();
650  } // *MENU*
651  virtual void SetParName(Int_t ipar, const char *name);
652  virtual void SetParNames(const char *name0 = "p0", const char *name1 = "p1", const char *name2 = "p2",
653  const char *name3 = "p3", const char *name4 = "p4", const char *name5 = "p5",
654  const char *name6 = "p6", const char *name7 = "p7", const char *name8 = "p8",
655  const char *name9 = "p9", const char *name10 = "p10"); // *MENU*
656  virtual void SetParError(Int_t ipar, Double_t error);
657  virtual void SetParErrors(const Double_t *errors);
658  virtual void SetParLimits(Int_t ipar, Double_t parmin, Double_t parmax);
659  virtual void SetParent(TObject *p = 0)
660  {
661  fParent = p;
662  }
663  virtual void SetRange(Double_t xmin, Double_t xmax); // *MENU*
664  virtual void SetRange(Double_t xmin, Double_t ymin, Double_t xmax, Double_t ymax);
665  virtual void SetRange(Double_t xmin, Double_t ymin, Double_t zmin, Double_t xmax, Double_t ymax, Double_t zmax);
666  virtual void SetSavedPoint(Int_t point, Double_t value);
667  virtual void SetTitle(const char *title = ""); // *MENU*
668  virtual void SetVectorized(Bool_t vectorized)
669  {
670  if (fType == EFType::kFormula && fFormula)
671  fFormula->SetVectorized(vectorized);
672  else
673  Warning("SetVectorized", "Can only set vectorized flag on formula-based TF1");
674  }
675  virtual void Update();
676 
677  static TF1 *GetCurrent();
678  static void AbsValue(Bool_t reject = kTRUE);
679  static void RejectPoint(Bool_t reject = kTRUE);
680  static Bool_t RejectedPoint();
681  static void SetCurrent(TF1 *f1);
682 
683  //Moments
684  virtual Double_t Moment(Double_t n, Double_t a, Double_t b, const Double_t *params = 0, Double_t epsilon = 0.000001);
685  virtual Double_t CentralMoment(Double_t n, Double_t a, Double_t b, const Double_t *params = 0, Double_t epsilon = 0.000001);
686  virtual Double_t Mean(Double_t a, Double_t b, const Double_t *params = 0, Double_t epsilon = 0.000001)
687  {
688  return Moment(1, a, b, params, epsilon);
689  }
690  virtual Double_t Variance(Double_t a, Double_t b, const Double_t *params = 0, Double_t epsilon = 0.000001)
691  {
692  return CentralMoment(2, a, b, params, epsilon);
693  }
694 
695  //some useful static utility functions to compute sampling points for Integral
696  //static void CalcGaussLegendreSamplingPoints(TGraph *g, Double_t eps=3.0e-11);
697  //static TGraph *CalcGaussLegendreSamplingPoints(Int_t num=21, Double_t eps=3.0e-11);
698  static void CalcGaussLegendreSamplingPoints(Int_t num, Double_t *x, Double_t *w, Double_t eps = 3.0e-11);
699 
700 private:
701  template <class T>
702  T EvalParTempl(const T *data, const Double_t *params = 0);
703 
704 #ifdef R__HAS_VECCORE
705  inline double EvalParVec(const Double_t *data, const Double_t *params);
706 #endif
707 
708  ClassDef(TF1, 10) // The Parametric 1-D function
709 };
710 
711 namespace ROOT {
712  namespace Internal {
713 
714  template<class Func>
715  void TF1Builder<Func>::Build(TF1 *f, Func func)
716  {
717  using Fnc_t = typename ROOT::Internal::GetFunctorType<decltype(ROOT::Internal::GetTheRightOp(&Func::operator()))>::type;
718  f->fType = std::is_same<Fnc_t, double>::value? TF1::EFType::kTemplScalar : TF1::EFType::kTemplVec;
719  f->fFunctor = new TF1::TF1FunctorPointerImpl<Fnc_t>(ROOT::Math::ParamFunctorTempl<Fnc_t>(func));
720  f->fParams = new TF1Parameters(f->fNpar);
721  }
722 
723  template<class Func>
724  void TF1Builder<Func *>::Build(TF1 *f, Func *func)
725  {
726  using Fnc_t = typename ROOT::Internal::GetFunctorType<decltype(ROOT::Internal::GetTheRightOp(&Func::operator()))>::type;
727  f->fType = std::is_same<Fnc_t, double>::value? TF1::EFType::kTemplScalar : TF1::EFType::kTemplVec;
728  f->fFunctor = new TF1::TF1FunctorPointerImpl<Fnc_t>(ROOT::Math::ParamFunctorTempl<Fnc_t>(func));
729  f->fParams = new TF1Parameters(f->fNpar);
730  }
731 
732  /// TF1 building from a string
733  /// used to build a TFormula based on a lambda function
734  template<>
735  struct TF1Builder<const char *> {
736  static void Build(TF1 *f, const char *formula)
737  {
738  f->fType = TF1::EFType::kFormula;
739  f->fFormula = new TFormula("tf1lambda", formula, f->fNdim, f->fNpar, false);
740  TString formulaExpression(formula);
741  Ssiz_t first = formulaExpression.Index("return") + 7;
742  Ssiz_t last = formulaExpression.Last(';');
743  TString title = formulaExpression(first, last - first);
744  f->SetTitle(title);
745  }
746  };
747  }
748 }
749 
750 inline Double_t TF1::operator()(Double_t x, Double_t y, Double_t z, Double_t t) const
751 {
752  return Eval(x, y, z, t);
753 }
754 
755 template <class T>
756 inline T TF1::operator()(const T *x, const Double_t *params)
757 {
758  return EvalPar(x, params);
759 }
760 
761 ////////////////////////////////////////////////////////////////////////////////
762 /// EvalPar for vectorized
763 template <class T>
764 T TF1::EvalPar(const T *x, const Double_t *params)
765 {
766  if (fType == EFType::kTemplVec || fType == EFType::kTemplScalar) {
767  return EvalParTempl(x, params);
768  } else if (fType == EFType::kFormula) {
769  return fFormula->EvalPar(x, params);
770  } else
771  return TF1::EvalPar((double *)x, params);
772 }
773 
774 ////////////////////////////////////////////////////////////////////////////////
775 /// Eval for vectorized functions
776 // template <class T>
777 // T TF1::Eval(T x, T y, T z, T t) const
778 // {
779 // if (fType == EFType::kFormula)
780 // return fFormula->Eval(x, y, z, t);
781 
782 // T xx[] = {x, y, z, t};
783 // Double_t *pp = (Double_t *)fParams->GetParameters();
784 // return ((TF1 *)this)->EvalPar(xx, pp);
785 // }
786 
787 // Internal to TF1. Evaluates Templated interfaces
788 template <class T>
789 inline T TF1::EvalParTempl(const T *data, const Double_t *params)
790 {
791  assert(fType == EFType::kTemplScalar || fType == EFType::kTemplVec);
792  if (!params) params = (Double_t *)fParams->GetParameters();
793  if (fFunctor)
794  return ((TF1FunctorPointerImpl<T> *)fFunctor)->fImpl(data, params);
795 
796  // this should throw an error
797  // we nned to implement a vectorized GetSave(x)
798  return TMath::SignalingNaN();
799 }
800 
801 #ifdef R__HAS_VECCORE
802 // Internal to TF1. Evaluates Vectorized TF1 on data of type Double_v
803 inline double TF1::EvalParVec(const Double_t *data, const Double_t *params)
804 {
805  assert(fType == EFType::kTemplVec);
806  std::vector<ROOT::Double_v> d(fNdim);
807  ROOT::Double_v res;
808 
809  for(auto i=0; i<fNdim; i++) {
810  d[i] = ROOT::Double_v(data[i]);
811  }
812 
813  if (fFunctor) {
814  res = ((TF1FunctorPointerImpl<ROOT::Double_v> *) fFunctor)->fImpl(d.data(), params);
815  } else {
816  // res = GetSave(x);
817  return TMath::SignalingNaN();
818  }
819  return vecCore::Get<ROOT::Double_v>(res, 0);
820 }
821 #endif
822 
823 inline void TF1::SetRange(Double_t xmin, Double_t, Double_t xmax, Double_t)
824 {
825  TF1::SetRange(xmin, xmax);
826 }
827 inline void TF1::SetRange(Double_t xmin, Double_t, Double_t, Double_t xmax, Double_t, Double_t)
828 {
829  TF1::SetRange(xmin, xmax);
830 }
831 
832 template <typename Func>
833 void TF1::SetFunction(Func f)
834 {
835  // set function from a generic C++ callable object
836  fType = EFType::kPtrScalarFreeFcn;
837  fFunctor = new TF1::TF1FunctorPointerImpl<double>(ROOT::Math::ParamFunctor(f));
838 }
839 template <class PtrObj, typename MemFn>
840 void TF1::SetFunction(PtrObj &p, MemFn memFn)
841 {
842  // set from a pointer to a member function
843  fType = EFType::kPtrScalarFreeFcn;
844  fFunctor = new TF1::TF1FunctorPointerImpl<double>(ROOT::Math::ParamFunctor(p, memFn));
845 }
846 
847 template <class T>
848 inline T TF1::GradientPar(Int_t ipar, const T *x, Double_t eps)
849 {
850  if (fType == EFType::kTemplVec || fType == EFType::kTemplScalar) {
851  return GradientParTempl<T>(ipar, x, eps);
852  } else
853  return GradientParTempl<Double_t>(ipar, (const Double_t *)x, eps);
854 }
855 
856 template <class T>
857 inline T TF1::GradientParTempl(Int_t ipar, const T *x, Double_t eps)
858 {
859  if (GetNpar() == 0)
860  return 0;
861 
862  if (eps < 1e-10 || eps > 1) {
863  Warning("Derivative", "parameter esp=%g out of allowed range[1e-10,1], reset to 0.01", eps);
864  eps = 0.01;
865  }
866  Double_t h;
867  TF1 *func = (TF1 *)this;
868  Double_t *parameters = GetParameters();
869 
870  // Copy parameters for thread safety
871  std::vector<Double_t> parametersCopy(parameters, parameters + GetNpar());
872  parameters = parametersCopy.data();
873 
874  Double_t al, bl, h2;
875  T f1, f2, g1, g2, d0, d2;
876 
877  ((TF1 *)this)->GetParLimits(ipar, al, bl);
878  if (al * bl != 0 && al >= bl) {
879  // this parameter is fixed
880  return 0;
881  }
882 
883  // check if error has been computer (is not zero)
884  if (func->GetParError(ipar) != 0)
885  h = eps * func->GetParError(ipar);
886  else
887  h = eps;
888 
889  // save original parameters
890  Double_t par0 = parameters[ipar];
891 
892  parameters[ipar] = par0 + h;
893  f1 = func->EvalPar(x, parameters);
894  parameters[ipar] = par0 - h;
895  f2 = func->EvalPar(x, parameters);
896  parameters[ipar] = par0 + h / 2;
897  g1 = func->EvalPar(x, parameters);
898  parameters[ipar] = par0 - h / 2;
899  g2 = func->EvalPar(x, parameters);
900 
901  // compute the central differences
902  h2 = 1 / (2. * h);
903  d0 = f1 - f2;
904  d2 = 2 * (g1 - g2);
905 
906  T grad = h2 * (4 * d2 - d0) / 3.;
907 
908  // restore original value
909  parameters[ipar] = par0;
910 
911  return grad;
912 }
913 
914 template <class T>
915 inline void TF1::GradientPar(const T *x, T *grad, Double_t eps)
916 {
917  if (fType == EFType::kTemplVec || fType == EFType::kTemplScalar) {
918  GradientParTempl<T>(x, grad, eps);
919  } else
920  GradientParTempl<Double_t>((const Double_t *)x, (Double_t *)grad, eps);
921 }
922 
923 template <class T>
924 inline void TF1::GradientParTempl(const T *x, T *grad, Double_t eps)
925 {
926  if (eps < 1e-10 || eps > 1) {
927  Warning("Derivative", "parameter esp=%g out of allowed range[1e-10,1], reset to 0.01", eps);
928  eps = 0.01;
929  }
930 
931  for (Int_t ipar = 0; ipar < GetNpar(); ipar++) {
932  grad[ipar] = GradientParTempl<T>(ipar, x, eps);
933  }
934 }
935 
936 #endif