Logo ROOT   6.08/07
Reference Guide
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 //////////////////////////////////////////////////////////////////////////
19 // //
20 // TF1 //
21 // //
22 // The Parametric 1-D function //
23 // //
24 //////////////////////////////////////////////////////////////////////////
25 
26 #include "RConfigure.h"
27 
28 #ifndef ROOT_TFormula
29 #include "TFormula.h"
30 #endif
31 #ifndef ROOT_TAttLine
32 #include "TAttLine.h"
33 #endif
34 #ifndef ROOT_TAttFill
35 #include "TAttFill.h"
36 #endif
37 #ifndef ROOT_TAttMarker
38 #include "TAttMarker.h"
39 #endif
40 
41 #ifndef ROOT_Math_ParamFunctor
42 #include "Math/ParamFunctor.h"
43 #endif
44 
45 class TF1;
46 class TH1;
47 class TAxis;
48 class TMethodCall;
49 
50 namespace ROOT {
51  namespace Fit {
52  class FitResult;
53  }
54 }
55 
57 public:
58  TF1Parameters() {} // needed for the I/O
60  fParameters(std::vector<Double_t>(npar) ),
61  fParNames(std::vector<std::string>(npar) )
62  {
63  for (int i = 0; i < npar; ++i ){
64  fParNames[i] = std::string(TString::Format("p%d",i).Data() );
65  }
66  }
67  // copy constructor
69  fParameters(rhs.fParameters),
70  fParNames(rhs.fParNames)
71  {}
72  // assignment
74  if (&rhs == this) return *this;
75  fParameters = rhs.fParameters;
76  fParNames = rhs.fParNames;
77  return *this;
78  }
79  virtual ~TF1Parameters() {}
80 
81  // getter methods
82  Double_t GetParameter(Int_t iparam) const {
83  return (CheckIndex(iparam)) ? fParameters[iparam] : 0;
84  }
85  Double_t GetParameter(const char * name) const {
86  return GetParameter(GetParNumber(name) );
87  }
88  const Double_t *GetParameters() const {
89  return fParameters.data();
90  }
91  const std::vector<double> & ParamsVec() const { return fParameters; }
92 
93  Int_t GetParNumber(const char * name) const;
94 
95  const char *GetParName(Int_t iparam) const {
96  return (CheckIndex(iparam)) ? fParNames[iparam].c_str() : "";
97  }
98 
99 
100  // setter methods
101  void SetParameter(Int_t iparam, Double_t value) {
102  if (!CheckIndex(iparam)) return;
103  fParameters[iparam] = value;
104  }
105  void SetParameters(const Double_t *params) {
106  std::copy(params, params + fParameters.size(), fParameters.begin() );
107  }
109  Double_t p5=0,Double_t p6=0,Double_t p7=0,Double_t p8=0,
110  Double_t p9=0,Double_t p10=0);
111 
112  void SetParameter(const char *name, Double_t value) {
113  SetParameter(GetParNumber(name),value);
114  }
115  void SetParName(Int_t iparam, const char * name) {
116  if (!CheckIndex(iparam)) return;
117  fParNames[iparam] = std::string(name);
118  }
119  void SetParNames(const char *name0="p0",const char *name1="p1",const char *name2="p2",
120  const char *name3="p3",const char*name4="p4", const char *name5="p5",
121  const char *name6="p6",const char *name7="p7",const char *name8="p8",
122  const char *name9="p9",const char *name10="p10");
123 
124 
125 
126  ClassDef(TF1Parameters,1) // The Parameters of a parameteric function
127 private:
128 
129  bool CheckIndex(Int_t i) const {
130  return (i >= 0 && i < int(fParameters.size()) );
131  }
132 
133  std::vector<Double_t> fParameters; // parameter values
134  std::vector<std::string> fParNames; // parameter names
135 };
136 
137 namespace ROOT {
138  namespace Internal {
139  /// Internal class used by TF1 for defining
140  /// template specialization for different TF1 constructors
141  template<class Func>
142  struct TF1Builder {
143  static void Build(TF1 * f, Func func);
144  };
145  }
146 }
147 
148 
149 class TF1 : public TNamed, public TAttLine, public TAttFill, public TAttMarker {
150 
151  template<class Func>
153 
154 public:
155  // Add to list behavior
156  enum class EAddToList {
157  kDefault,
158  kAdd,
159  kNo
160  };
161 
162 protected:
163  Double_t fXmin; //Lower bounds for the range
164  Double_t fXmax; //Upper bounds for the range
165  Int_t fNpar; //Number of parameters
166  Int_t fNdim; //Function dimension
167  Int_t fNpx; //Number of points used for the graphical representation
168  Int_t fType; //(=0 for standard functions, 1 if pointer to function)
169  Int_t fNpfits; //Number of points used in the fit
170  Int_t fNDF; //Number of degrees of freedom in the fit
171  Double_t fChisquare; //Function fit chisquare
172  Double_t fMinimum; //Minimum value for plotting
173  Double_t fMaximum; //Maximum value for plotting
174  std::vector<Double_t> fParErrors; //Array of errors of the fNpar parameters
175  std::vector<Double_t> fParMin; //Array of lower limits of the fNpar parameters
176  std::vector<Double_t> fParMax; //Array of upper limits of the fNpar parameters
177  std::vector<Double_t> fSave; //Array of fNsave function values
178  std::vector<Double_t> fIntegral; //!Integral of function binned on fNpx bins
179  std::vector<Double_t> fAlpha; //!Array alpha. for each bin in x the deconvolution r of fIntegral
180  std::vector<Double_t> fBeta; //!Array beta. is approximated by x = alpha +beta*r *gamma*r**2
181  std::vector<Double_t> fGamma; //!Array gamma.
182  TObject *fParent; //!Parent object hooking this function (if one)
183  TH1 *fHistogram; //!Pointer to histogram used for visualisation
184  TMethodCall *fMethodCall; //!Pointer to MethodCall in case of interpreted function
185  Bool_t fNormalized; //Normalization option (false by default)
186  Double_t fNormIntegral;//Integral of the function before being normalized
187  ROOT::Math::ParamFunctor fFunctor; //! Functor object to wrap any C++ callable object
188  TFormula *fFormula; //Pointer to TFormula in case when user define formula
189  TF1Parameters *fParams; //Pointer to Function parameters object (exusts only for not-formula functions)
190 
191  static std::atomic<Bool_t> fgAbsValue; //use absolute value of function when computing integral
192  static Bool_t fgRejectPoint; //True if point must be rejected in a fit
193  static std::atomic<Bool_t> fgAddToGlobList; //True if we want to register the function in the global list
194  static TF1 *fgCurrent; //pointer to current function being processed
195 
196 
197  //void CreateFromFunctor(const char *name, Int_t npar, Int_t ndim = 1);
198  void DoInitialize(EAddToList addToGlobList);
199 
200  void IntegrateForNormalization();
201 
202  virtual Double_t GetMinMaxNDim(Double_t * x , Bool_t findmax, Double_t epsilon = 0, Int_t maxiter = 0) const;
203  virtual void GetRange(Double_t * xmin, Double_t * xmax) const;
204  virtual TH1 *DoCreateHistogram(Double_t xmin, Double_t xmax, Bool_t recreate = kFALSE);
205 
206  enum {
207  kNotGlobal = BIT(10), // don't register in global list of functions
208  };
209 
210 public:
211  // TF1 status bits
212  enum {
213  kNotDraw = BIT(9) // don't draw the function when in a TH1
214  };
215 
216  TF1();
217  TF1(const char *name, const char *formula, Double_t xmin=0, Double_t xmax = 1, EAddToList addToGlobList = EAddToList::kDefault);
218  TF1(const char *name, Double_t xmin, Double_t xmax, Int_t npar,Int_t ndim = 1, EAddToList addToGlobList = EAddToList::kDefault);
219  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);
220  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);
221 
222  // Constructors using functors (compiled mode only)
223  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);
224 
225  // Template constructors from any C++ callable object, defining the operator() (double * , double *)
226  // and returning a double.
227  // The class name is not needed when using compile code, while it is required when using
228  // interpreted code via the specialized constructor with void *.
229  // An instance of the C++ function class or its pointer can both be used. The former is reccomended when using
230  // C++ compiled code, but if CINT compatibility is needed, then a pointer to the function class must be used.
231  // xmin and xmax specify the plotting range, npar is the number of parameters.
232  // See the tutorial math/exampleFunctor.C for an example of using this constructor
233  template <typename Func>
234  TF1(const char *name, Func f, Double_t xmin, Double_t xmax, Int_t npar,Int_t ndim = 1, EAddToList addToGlobList = EAddToList::kDefault );
235 
236  // backward compatible interface
237  template <typename Func>
238  TF1(const char *name, Func f, Double_t xmin, Double_t xmax, Int_t npar, const char *, EAddToList addToGlobList = EAddToList::kDefault ) :
239  TNamed(name,name), TAttLine(), TAttFill(), TAttMarker(),
240  fXmin(xmin), fXmax(xmax),
241  fNpar(npar), fNdim(1),
242  fNpx(100), fType(1),
243  fNpfits(0), fNDF(0), fChisquare(0),
244  fMinimum(-1111), fMaximum(-1111),
245  fParErrors(std::vector<Double_t>(npar)),
246  fParMin(std::vector<Double_t>(npar)),
247  fParMax(std::vector<Double_t>(npar)),
248  fParent(0), fHistogram(0),
249  fMethodCall(0),
250  fNormalized(false), fNormIntegral(0),
251  fFunctor(ROOT::Math::ParamFunctor(f)),
252  fFormula(0),
253  fParams(new TF1Parameters(npar) )
254  {
255  DoInitialize(addToGlobList);
256  }
257 
258 
259  // Template constructors from a pointer to any C++ class of type PtrObj with a specific member function of type
260  // MemFn.
261  // The member function must have the signature of (double * , double *) and returning a double.
262  // The class name and the method name are not needed when using compile code
263  // (the member function pointer is used in this case), while they are required when using interpreted
264  // code via the specialized constructor with void *.
265  // xmin and xmax specify the plotting range, npar is the number of parameters.
266  // See the tutorial math/exampleFunctor.C for an example of using this constructor
267  template <class PtrObj, typename MemFn>
268  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) :
269  TNamed(name,name), TAttLine(), TAttFill(), TAttMarker(),
270  fXmin(xmin), fXmax(xmax),
271  fNpar(npar), fNdim(ndim),
272  fNpx(100), fType(1),
273  fNpfits(0), fNDF(0), fChisquare(0),
274  fMinimum(-1111), fMaximum(-1111),
275  fParErrors(std::vector<Double_t>(npar)),
276  fParMin(std::vector<Double_t>(npar)),
277  fParMax(std::vector<Double_t>(npar)),
278  fParent(0), fHistogram(0),
279  fMethodCall(0),
280  fNormalized(false), fNormIntegral(0),
281  fFunctor ( ROOT::Math::ParamFunctor(p,memFn) ),
282  fFormula(0),
283  fParams(new TF1Parameters(npar) )
284  {
285  DoInitialize(addToGlobList);
286  }
287  // backward compatible interface
288  template <class PtrObj, typename MemFn>
289  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) :
290  TNamed(name,name), TAttLine(), TAttFill(), TAttMarker(),
291  fXmin(xmin), fXmax(xmax),
292  fNpar(npar), fNdim(1),
293  fNpx(100), fType(1),
294  fNpfits(0), fNDF(0), fChisquare(0),
295  fMinimum(-1111), fMaximum(-1111),
296  fParErrors(std::vector<Double_t>(npar)),
297  fParMin(std::vector<Double_t>(npar)),
298  fParMax(std::vector<Double_t>(npar)),
299  fParent(0), fHistogram(0),
300  fMethodCall(0),
301  fNormalized(false), fNormIntegral(0),
302  fFunctor ( ROOT::Math::ParamFunctor(p,memFn) ),
303  fFormula(0),
304  fParams(new TF1Parameters(npar) )
305  {
306  DoInitialize(addToGlobList);
307  }
308 
309  TF1(const TF1 &f1);
310  TF1& operator=(const TF1 &rhs);
311  virtual ~TF1();
312  virtual void AddParameter(const TString &name, Double_t value) { if (fFormula) fFormula->AddParameter(name,value); }
313  // virtual void AddParameters(const pair<TString,Double_t> *pairs, Int_t size) { fFormula->AddParameters(pairs,size); }
314  // virtual void AddVariable(const TString &name, Double_t value = 0) { if (fFormula) fFormula->AddVariable(name,value); }
315  // virtual void AddVariables(const TString *vars, Int_t size) { if (fFormula) fFormula->AddVariables(vars,size); }
316  virtual Bool_t AddToGlobalList(Bool_t on = kTRUE);
317  static Bool_t DefaultAddToGlobalList(Bool_t on = kTRUE);
318  virtual void Browse(TBrowser *b);
319  virtual void Copy(TObject &f1) const;
320  virtual Double_t Derivative (Double_t x, Double_t *params=0, Double_t epsilon=0.001) const;
321  virtual Double_t Derivative2(Double_t x, Double_t *params=0, Double_t epsilon=0.001) const;
322  virtual Double_t Derivative3(Double_t x, Double_t *params=0, Double_t epsilon=0.001) const;
323  static Double_t DerivativeError();
324  virtual Int_t DistancetoPrimitive(Int_t px, Int_t py);
325  virtual void Draw(Option_t *option="");
326  virtual TF1 *DrawCopy(Option_t *option="") const;
327  virtual TObject *DrawDerivative(Option_t *option="al"); // *MENU*
328  virtual TObject *DrawIntegral(Option_t *option="al"); // *MENU*
329  virtual void DrawF1(Double_t xmin, Double_t xmax, Option_t *option="");
330  virtual Double_t Eval(Double_t x, Double_t y=0, Double_t z=0, Double_t t=0) const;
331  virtual Double_t EvalPar(const Double_t *x, const Double_t *params=0);
332  virtual Double_t operator()(Double_t x, Double_t y=0, Double_t z = 0, Double_t t = 0) const;
333  virtual Double_t operator()(const Double_t *x, const Double_t *params=0);
334  virtual void ExecuteEvent(Int_t event, Int_t px, Int_t py);
335  virtual void FixParameter(Int_t ipar, Double_t value);
336  Double_t GetChisquare() const {return fChisquare;}
337  virtual TH1 *GetHistogram() const;
338  virtual TH1 *CreateHistogram() { return DoCreateHistogram(fXmin, fXmax); }
339  virtual TFormula *GetFormula() { return fFormula;}
340  virtual const TFormula *GetFormula() const { return fFormula;}
341  virtual TString GetExpFormula(Option_t *option="") const { return (fFormula) ? fFormula->GetExpFormula(option) : ""; }
342  virtual const TObject *GetLinearPart(Int_t i) const { return (fFormula) ? fFormula->GetLinearPart(i) : nullptr;}
343  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;
344  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;
345  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;
346  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;
347  virtual Double_t GetMaximumStored() const {return fMaximum;}
348  virtual Double_t GetMinimumStored() const {return fMinimum;}
349  virtual Int_t GetNpar() const { return fNpar;}
350  virtual Int_t GetNdim() const { return fNdim;}
351  virtual Int_t GetNDF() const;
352  virtual Int_t GetNpx() const {return fNpx;}
353  TMethodCall *GetMethodCall() const {return fMethodCall;}
354  virtual Int_t GetNumber() const { return (fFormula) ? fFormula->GetNumber() : 0;}
355  virtual Int_t GetNumberFreeParameters() const;
356  virtual Int_t GetNumberFitPoints() const {return fNpfits;}
357  virtual char *GetObjectInfo(Int_t px, Int_t py) const;
358  TObject *GetParent() const {return fParent;}
359  virtual Double_t GetParameter(Int_t ipar) const {
360  return (fFormula) ? fFormula->GetParameter(ipar) : fParams->GetParameter(ipar);
361  }
362  virtual Double_t GetParameter(const TString &name) const {
363  return (fFormula) ? fFormula->GetParameter(name) : fParams->GetParameter(name);
364  }
365  virtual Double_t *GetParameters() const {
366  return (fFormula) ? fFormula->GetParameters() : const_cast<Double_t*>(fParams->GetParameters());
367  }
368  virtual void GetParameters(Double_t *params) { if (fFormula) fFormula->GetParameters(params);
369  else std::copy(fParams->ParamsVec().begin(), fParams->ParamsVec().end(), params); }
370  virtual const char *GetParName(Int_t ipar) const {
371  return (fFormula) ? fFormula->GetParName(ipar) : fParams->GetParName(ipar);
372  }
373  virtual Int_t GetParNumber(const char* name) const {
374  return (fFormula) ? fFormula->GetParNumber(name) : fParams->GetParNumber(name);
375  }
376  virtual Double_t GetParError(Int_t ipar) const;
377  virtual const Double_t *GetParErrors() const {return fParErrors.data();}
378  virtual void GetParLimits(Int_t ipar, Double_t &parmin, Double_t &parmax) const;
379  virtual Double_t GetProb() const;
380  virtual Int_t GetQuantiles(Int_t nprobSum, Double_t *q, const Double_t *probSum);
381  virtual Double_t GetRandom();
382  virtual Double_t GetRandom(Double_t xmin, Double_t xmax);
383  virtual void GetRange(Double_t &xmin, Double_t &xmax) const;
384  virtual void GetRange(Double_t &xmin, Double_t &ymin, Double_t &xmax, Double_t &ymax) const;
385  virtual void GetRange(Double_t &xmin, Double_t &ymin, Double_t &zmin, Double_t &xmax, Double_t &ymax, Double_t &zmax) const;
386  virtual Double_t GetSave(const Double_t *x);
387  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;
388  virtual Double_t GetXmin() const {return fXmin;}
389  virtual Double_t GetXmax() const {return fXmax;}
390  TAxis *GetXaxis() const ;
391  TAxis *GetYaxis() const ;
392  TAxis *GetZaxis() const ;
393  virtual Double_t GetVariable(const TString &name) { return (fFormula) ? fFormula->GetVariable(name) : 0;}
394  virtual Double_t GradientPar(Int_t ipar, const Double_t *x, Double_t eps=0.01);
395  virtual void GradientPar(const Double_t *x, Double_t *grad, Double_t eps=0.01);
396  virtual void InitArgs(const Double_t *x, const Double_t *params);
397  static void InitStandardFunctions();
398  virtual Double_t Integral(Double_t a, Double_t b, Double_t epsrel=1.e-12);
399  virtual Double_t IntegralOneDim(Double_t a, Double_t b, Double_t epsrel, Double_t epsabs, Double_t &err);
400  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);
401  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);
402  // virtual Double_t IntegralFast(const TGraph *g, Double_t a, Double_t b, Double_t *params=0);
403  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);
404  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);
405  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) {
406  return IntegralMultiple(n,a,b,maxpts, epsrel, epsrel, relerr, nfnevl, ifail);
407  }
408  virtual Double_t IntegralMultiple(Int_t n, const Double_t *a, const Double_t *b, Double_t epsrel, Double_t &relerr);
409  virtual Bool_t IsEvalNormalized() const { return fNormalized; }
410  /// return kTRUE if the point is inside the function range
411  virtual Bool_t IsInside(const Double_t *x) const { return !( ( x[0] < fXmin) || ( x[0] > fXmax ) ); }
412  virtual Bool_t IsLinear() const { return (fFormula) ? fFormula->IsLinear() : false;}
413  virtual Bool_t IsValid() const;
414  virtual void Print(Option_t *option="") const;
415  virtual void Paint(Option_t *option="");
416  virtual void ReleaseParameter(Int_t ipar);
417  virtual void Save(Double_t xmin, Double_t xmax, Double_t ymin, Double_t ymax, Double_t zmin, Double_t zmax);
418  virtual void SavePrimitive(std::ostream &out, Option_t *option = "");
419  virtual void SetChisquare(Double_t chi2) {fChisquare = chi2;}
420  virtual void SetFitResult(const ROOT::Fit::FitResult & result, const Int_t * indpar = 0);
421  template <class PtrObj, typename MemFn>
422  void SetFunction( PtrObj& p, MemFn memFn );
423  template <typename Func>
424  void SetFunction( Func f );
425  virtual void SetMaximum(Double_t maximum=-1111); // *MENU*
426  virtual void SetMinimum(Double_t minimum=-1111); // *MENU*
427  virtual void SetNDF(Int_t ndf);
428  virtual void SetNumberFitPoints(Int_t npfits) {fNpfits = npfits;}
429  virtual void SetNormalized(Bool_t flag) { fNormalized = flag; Update(); }
430  virtual void SetNpx(Int_t npx=100); // *MENU*
431  virtual void SetParameter(Int_t param, Double_t value) {
432  (fFormula) ? fFormula->SetParameter(param,value) : fParams->SetParameter(param,value);
433  Update();
434  }
435  virtual void SetParameter(const TString &name, Double_t value) {
436  (fFormula) ? fFormula->SetParameter(name,value) : fParams->SetParameter(name,value);
437  Update();
438  }
439  virtual void SetParameters(const Double_t *params) {
440  (fFormula) ? fFormula->SetParameters(params) : fParams->SetParameters(params);
441  Update();
442  }
444  Double_t p5=0,Double_t p6=0,Double_t p7=0,Double_t p8=0,
445  Double_t p9=0,Double_t p10=0) {
446  if (fFormula) fFormula->SetParameters(p0,p1,p2,p3,p4,p5,p6,p7,p8,p9,p10);
447  else fParams->SetParameters(p0,p1,p2,p3,p4,p5,p6,p7,p8,p9,p10);
448  Update();
449  } // *MENU*
450  virtual void SetParName(Int_t ipar, const char *name);
451  virtual void SetParNames(const char *name0="p0",const char *name1="p1",const char *name2="p2",
452  const char *name3="p3",const char *name4="p4", const char *name5="p5",
453  const char *name6="p6",const char *name7="p7",const char *name8="p8",
454  const char *name9="p9",const char *name10="p10"); // *MENU*
455  virtual void SetParError(Int_t ipar, Double_t error);
456  virtual void SetParErrors(const Double_t *errors);
457  virtual void SetParLimits(Int_t ipar, Double_t parmin, Double_t parmax);
458  virtual void SetParent(TObject *p=0) {fParent = p;}
459  virtual void SetRange(Double_t xmin, Double_t xmax); // *MENU*
460  virtual void SetRange(Double_t xmin, Double_t ymin, Double_t xmax, Double_t ymax);
461  virtual void SetRange(Double_t xmin, Double_t ymin, Double_t zmin, Double_t xmax, Double_t ymax, Double_t zmax);
462  virtual void SetSavedPoint(Int_t point, Double_t value);
463  virtual void SetTitle(const char *title=""); // *MENU*
464  virtual void Update();
465 
466  static TF1 *GetCurrent();
467  static void AbsValue(Bool_t reject=kTRUE);
468  static void RejectPoint(Bool_t reject=kTRUE);
469  static Bool_t RejectedPoint();
470  static void SetCurrent(TF1 *f1);
471 
472  //Moments
473  virtual Double_t Moment(Double_t n, Double_t a, Double_t b, const Double_t *params=0, Double_t epsilon=0.000001);
474  virtual Double_t CentralMoment(Double_t n, Double_t a, Double_t b, const Double_t *params=0, Double_t epsilon=0.000001);
475  virtual Double_t Mean(Double_t a, Double_t b, const Double_t *params=0, Double_t epsilon=0.000001) {return Moment(1,a,b,params,epsilon);}
476  virtual Double_t Variance(Double_t a, Double_t b, const Double_t *params=0, Double_t epsilon=0.000001) {return CentralMoment(2,a,b,params,epsilon);}
477 
478  //some useful static utility functions to compute sampling points for Integral
479  //static void CalcGaussLegendreSamplingPoints(TGraph *g, Double_t eps=3.0e-11);
480  //static TGraph *CalcGaussLegendreSamplingPoints(Int_t num=21, Double_t eps=3.0e-11);
481  static void CalcGaussLegendreSamplingPoints(Int_t num, Double_t *x, Double_t *w, Double_t eps=3.0e-11);
482 
483  ClassDef(TF1,9) //The Parametric 1-D function
484 };
485 
486 ///ctor implementation
487 template <typename Func>
488 TF1::TF1(const char *name, Func f, Double_t xmin, Double_t xmax, Int_t npar,Int_t ndim, EAddToList addToGlobList) :
489  TNamed(name,name), TAttLine(), TAttFill(), TAttMarker(),
490  fXmin(xmin), fXmax(xmax),
491  fNpar(npar), fNdim(ndim),
492  fNpx(100), fType(1),
493  fNpfits(0), fNDF(0), fChisquare(0),
494  fMinimum(-1111), fMaximum(-1111),
495  fParErrors(std::vector<Double_t>(npar)),
496  fParMin(std::vector<Double_t>(npar)),
497  fParMax(std::vector<Double_t>(npar)),
498  fParent(0), fHistogram(0),
499  fMethodCall(0),
500  fNormalized(false), fNormIntegral(0),
501  //fFunctor(ROOT::Math::ParamFunctor(f)),
502  fFormula(0),
503  fParams(0)
504 {
506  DoInitialize(addToGlobList);
507 }
508 
509 namespace ROOT {
510  namespace Internal {
511 
512  template<class Func>
514  f->fType = 1;
516  f->fParams = new TF1Parameters(f->fNpar);
517  }
518  /// TF1 building from a string
519  /// used to build a TFormula based on a lambda function
520  template<>
521  struct TF1Builder<const char *> {
522  static void Build(TF1 * f, const char * formula) {
523  f->fType = 0;
524  f->fFormula = new TFormula("tf1lambda",formula,f->fNdim,f->fNpar,false);
525  TString formulaExpression(formula);
526  Ssiz_t first = formulaExpression.Index("return")+7;
527  Ssiz_t last = formulaExpression.Last(';');
528  TString title = formulaExpression( first, last-first);
529  f->SetTitle(title);
530  }
531  };
532  }
533 }
534 
535 
536 
537 
539  { return Eval(x,y,z,t); }
540 inline Double_t TF1::operator()(const Double_t *x, const Double_t *params)
541  {
542 
543  if (fMethodCall) InitArgs(x,params);
544  return EvalPar(x,params);
545  }
546 
547 
549  { TF1::SetRange(xmin, xmax); }
551  { TF1::SetRange(xmin, xmax); }
552 
553 template <typename Func>
555  // set function from a generic C++ callable object
556  fType = 1;
557  fFunctor = ROOT::Math::ParamFunctor(f);
558 }
559 template <class PtrObj, typename MemFn>
560 void TF1::SetFunction( PtrObj& p, MemFn memFn ) {
561  // set from a pointer to a member function
562  fType = 1;
563  fFunctor = ROOT::Math::ParamFunctor(p,memFn);
564 }
565 
566 
567 
568 
569 #endif
ROOT::Math::ParamFunctor fFunctor
Definition: TF1.h:187
virtual TString GetExpFormula(Option_t *option="") const
Definition: TF1.h:341
Double_t fMaximum
Definition: TF1.h:173
virtual void SetParameters(Double_t p0, Double_t p1, Double_t p2=0, Double_t p3=0, Double_t p4=0, Double_t p5=0, Double_t p6=0, Double_t p7=0, Double_t p8=0, Double_t p9=0, Double_t p10=0)
Definition: TF1.h:443
TF1Parameters(Int_t npar)
Definition: TF1.h:59
virtual void AddParameter(const TString &name, Double_t value)
Definition: TF1.h:312
virtual void SetParameters(const Double_t *params)
Definition: TF1.h:439
Int_t fNpx
Definition: TF1.h:167
Bool_t fNormalized
Pointer to MethodCall in case of interpreted function.
Definition: TF1.h:185
float xmin
Definition: THbookFile.cxx:93
virtual Double_t operator()(Double_t x, Double_t y=0, Double_t z=0, Double_t t=0) const
Definition: TF1.h:538
virtual TFormula * GetFormula()
Definition: TF1.h:339
virtual Int_t GetNpx() const
Definition: TF1.h:352
static double p3(double t, double a, double b, double c, double d)
This namespace contains pre-defined functions to be used in conjuction with TExecutor::Map and TExecu...
Definition: StringConv.hxx:21
static std::atomic< Bool_t > fgAbsValue
Definition: TF1.h:191
const char Option_t
Definition: RtypesCore.h:62
float ymin
Definition: THbookFile.cxx:93
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)
Definition: TF1.h:289
TMethodCall * GetMethodCall() const
Definition: TF1.h:353
void SetParameters(const Double_t *params)
Definition: TFormula.cxx:2432
void DoInitialize(EAddToList addToGlobList)
Common initialization of the TF1.
Definition: TF1.cxx:611
#define BIT(n)
Definition: Rtypes.h:120
virtual Double_t GetMinimumStored() const
Definition: TF1.h:348
Int_t fNpar
Definition: TF1.h:165
TObject * fParent
Array gamma.
Definition: TF1.h:182
virtual void SetRange(Double_t xmin, Double_t xmax)
Initialize the upper and lower bounds to draw the function.
Definition: TF1.cxx:3240
Double_t fChisquare
Definition: TF1.h:171
Bool_t IsLinear() const
Definition: TFormula.h:187
Ssiz_t Index(const char *pat, Ssiz_t i=0, ECaseCompare cmp=kExact) const
Definition: TString.h:582
std::vector< std::string > fParNames
Definition: TF1.h:134
Basic string class.
Definition: TString.h:137
virtual void SetNumberFitPoints(Int_t npfits)
Definition: TF1.h:428
Double_t fNormIntegral
Definition: TF1.h:186
int Int_t
Definition: RtypesCore.h:41
void SetParameter(const char *name, Double_t value)
Definition: TF1.h:112
bool Bool_t
Definition: RtypesCore.h:59
TArc * a
Definition: textangle.C:12
const Bool_t kFALSE
Definition: Rtypes.h:92
void SetFunction(PtrObj &p, MemFn memFn)
Definition: TF1.h:560
Internal class used by TF1 for defining template specialization for different TF1 constructors...
Definition: TF1.h:142
STL namespace.
virtual void SetNormalized(Bool_t flag)
Definition: TF1.h:429
TF1 Parameters class.
Definition: TF1.h:56
Double_t GetParameter(Int_t iparam) const
Definition: TF1.h:82
Int_t GetNumber() const
Definition: TFormula.h:176
TF1(const char *name, Func f, Double_t xmin, Double_t xmax, Int_t npar, const char *, EAddToList addToGlobList=EAddToList::kDefault)
Definition: TF1.h:238
Marker Attributes class.
Definition: TAttMarker.h:24
virtual void SetParent(TObject *p=0)
Definition: TF1.h:458
TFormula * fFormula
Functor object to wrap any C++ callable object.
Definition: TF1.h:188
Fill Area Attributes class.
Definition: TAttFill.h:24
void SetParameters(TFitEditor::FuncParams_t &pars, TF1 *func)
Restore the parameters from pars into the function.
Definition: TFitEditor.cxx:287
Double_t x[n]
Definition: legend1.C:17
#define ClassDef(name, id)
Definition: Rtypes.h:254
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:2335
Int_t fNpfits
Definition: TF1.h:169
Int_t fNdim
Definition: TF1.h:166
The TNamed class is the base class for all named ROOT classes.
Definition: TNamed.h:33
const char * GetParName(Int_t iparam) const
Definition: TF1.h:95
virtual Double_t GetVariable(const TString &name)
Definition: TF1.h:393
virtual void GetParameters(Double_t *params)
Definition: TF1.h:368
bool CheckIndex(Int_t i) const
Definition: TF1.h:129
void SetParameter(const char *name, Double_t value)
Definition: TFormula.cxx:2352
static double p2(double t, double a, double b, double c)
TF1()
TF1 default constructor.
Definition: TF1.cxx:396
void SetParameters(const Double_t *params)
Definition: TF1.h:105
const TObject * GetLinearPart(Int_t i) const
Definition: TFormula.cxx:1999
std::vector< std::vector< double > > Data
virtual Double_t IntegralMultiple(Int_t n, const Double_t *a, const Double_t *b, Int_t, Int_t maxpts, Double_t epsrel, Double_t &relerr, Int_t &nfnevl, Int_t &ifail)
Definition: TF1.h:405
virtual Int_t GetNdim() const
Definition: TF1.h:350
Method or function calling interface.
Definition: TMethodCall.h:41
virtual Double_t GetMaximumStored() const
Definition: TF1.h:347
static void GetRange(const char *comments, Double_t &xmin, Double_t &xmax, Double_t &factor)
Parse comments to search for a range specifier of the style: [xmin,xmax] or [xmin,xmax,nbits] [0,1] [-10,100]; [-pi,pi], [-pi/2,pi/4],[-2pi,2*pi] [-10,100,16] [0,0,8] if nbits is not specified, or nbits <2 or nbits>32 it is set to 32 if (xmin==0 and xmax==0 and nbits <=16) the double word will be converted to a float and its mantissa truncated to nbits significative bits.
TF1Parameters & operator=(const TF1Parameters &rhs)
Definition: TF1.h:73
std::vector< Double_t > fIntegral
Definition: TF1.h:178
virtual const char * GetParName(Int_t ipar) const
Definition: TF1.h:370
static TF1 * fgCurrent
Definition: TF1.h:194
Double_t GetVariable(const char *name) const
Definition: TFormula.cxx:2142
Using a TBrowser one can browse all ROOT objects.
Definition: TBrowser.h:41
th1 Draw()
virtual void SetChisquare(Double_t chi2)
Definition: TF1.h:419
float ymax
Definition: THbookFile.cxx:93
static Bool_t fgRejectPoint
Definition: TF1.h:192
TF1Parameters()
Definition: TF1.h:58
std::vector< Double_t > fParameters
Definition: TF1.h:133
TH1 * fHistogram
Parent object hooking this function (if one)
Definition: TF1.h:183
virtual const Double_t * GetParErrors() const
Definition: TF1.h:377
Class to manage histogram axis.
Definition: TAxis.h:36
IParamFunction interface (abstract class) describing multi-dimensional parameteric functions It is a ...
virtual Int_t GetNumberFitPoints() const
Definition: TF1.h:356
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 F O R M U L A class.
Definition: TFormula.h:89
virtual const TObject * GetLinearPart(Int_t i) const
Definition: TF1.h:342
Int_t fType
Definition: TF1.h:168
virtual Bool_t IsLinear() const
Definition: TF1.h:412
Double_t E()
Definition: TMath.h:54
Double_t fMinimum
Definition: TF1.h:172
Double_t GetParameter(const char *name) const
Definition: TF1.h:85
virtual const TFormula * GetFormula() const
Definition: TF1.h:340
float xmax
Definition: THbookFile.cxx:93
static double p1(double t, double a, double b)
const std::vector< double > & ParamsVec() const
Definition: TF1.h:91
virtual Double_t GetParameter(const TString &name) const
Definition: TF1.h:362
Int_t GetParNumber(const char *name) const
Definition: TFormula.cxx:2282
virtual Double_t GetXmin() const
Definition: TF1.h:388
virtual void SetParameter(const TString &name, Double_t value)
Definition: TF1.h:435
REAL epsilon
Definition: triangle.c:617
void SetParameter(Int_t iparam, Double_t value)
Definition: TF1.h:101
virtual Double_t Mean(Double_t a, Double_t b, const Double_t *params=0, Double_t epsilon=0.000001)
Definition: TF1.h:475
virtual Bool_t IsEvalNormalized() const
Definition: TF1.h:409
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:3267
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)
Definition: TF1.h:268
int Ssiz_t
Definition: RtypesCore.h:63
Double_t GetChisquare() const
Definition: TF1.h:336
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:1196
class containg the result of the fit and all the related information (fitted parameter values...
Definition: FitResult.h:52
Int_t fNDF
Definition: TF1.h:170
void Copy(void *source, void *dest)
double f(double x)
void Print(std::ostream &os, const OptionType &opt)
double Double_t
Definition: RtypesCore.h:55
std::vector< Double_t > fSave
Definition: TF1.h:177
Ssiz_t Last(char c) const
Find last occurrence of a character c.
Definition: TString.cxx:865
std::vector< Double_t > fParErrors
Definition: TF1.h:174
Double_t fXmin
Definition: TF1.h:163
TFitResultPtr Fit(FitObject *h1, TF1 *f1, Foption_t &option, const ROOT::Math::MinimizerOptions &moption, const char *goption, ROOT::Fit::DataRange &range)
Definition: HFitImpl.cxx:134
void SetParName(Int_t iparam, const char *name)
Definition: TF1.h:115
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:3570
void AddParameter(const TString &name, Double_t value=0)
Definition: TFormula.h:161
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:181
double func(double *x, double *p)
Definition: stressTF1.cxx:213
const char * GetParName(Int_t ipar) const
Definition: TFormula.cxx:2319
virtual Double_t Variance(Double_t a, Double_t b, const Double_t *params=0, Double_t epsilon=0.000001)
Definition: TF1.h:476
virtual Double_t GetXmax() const
Definition: TF1.h:389
The TH1 histogram class.
Definition: TH1.h:80
you should not use this method at all Int_t Int_t Double_t Double_t Double_t e
Definition: TRolke.cxx:630
std::vector< Double_t > fAlpha
Integral of function binned on fNpx bins.
Definition: TF1.h:179
EAddToList
Definition: TF1.h:156
Namespace for new Math classes and functions.
virtual void InitArgs(const Double_t *x, const Double_t *params)
Initialize parameters addresses.
Definition: TF1.cxx:2234
std::vector< Double_t > fBeta
Array alpha. for each bin in x the deconvolution r of fIntegral.
Definition: TF1.h:180
TObject * GetParent() const
Definition: TF1.h:358
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
virtual Bool_t IsInside(const Double_t *x) const
return kTRUE if the point is inside the function range
Definition: TF1.h:411
TMethodCall * fMethodCall
Pointer to histogram used for visualisation.
Definition: TF1.h:184
virtual TH1 * CreateHistogram()
Definition: TF1.h:338
RooCmdArg Save(Bool_t flag=kTRUE)
virtual Int_t GetNpar() const
Definition: TF1.h:349
virtual ~TF1Parameters()
Definition: TF1.h:79
static void Build(TF1 *f, Func func)
Definition: TF1.h:513
TString GetExpFormula(Option_t *option="") const
return the expression formula If option = "P" replace the parameter names with their values If option...
Definition: TFormula.cxx:2641
virtual Double_t GetParameter(Int_t ipar) const
Definition: TF1.h:359
1-Dim function class
Definition: TF1.h:149
TF1 * f1
Definition: legend1.C:11
Param Functor class for Multidimensional functions.
Definition: ParamFunctor.h:209
static void Build(TF1 *f, const char *formula)
Definition: TF1.h:522
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
const Double_t * GetParameters() const
Definition: TF1.h:88
TF1Parameters * fParams
Definition: TF1.h:189
static std::atomic< Bool_t > fgAddToGlobList
Definition: TF1.h:193
virtual Double_t * GetParameters() const
Definition: TF1.h:365
Double_t GetParameter(const char *name) const
Definition: TFormula.cxx:2294
double result[121]
virtual void SetParameter(Int_t param, Double_t value)
Definition: TF1.h:431
Definition: first.py:1
std::vector< Double_t > fParMax
Definition: TF1.h:176
const Bool_t kTRUE
Definition: Rtypes.h:91
float * q
Definition: THbookFile.cxx:87
virtual Double_t EvalPar(const Double_t *x, const Double_t *params=0)
Evaluate function with given coordinates and parameters.
Definition: TF1.cxx:1225
virtual Int_t GetNumber() const
Definition: TF1.h:354
std::vector< Double_t > fParMin
Definition: TF1.h:175
const Int_t n
Definition: legend1.C:16
Line Attributes class.
Definition: TAttLine.h:24
Double_t fXmax
Definition: TF1.h:164
Double_t * GetParameters() const
Definition: TFormula.cxx:2334
char name[80]
Definition: TGX11.cxx:109
virtual Int_t GetParNumber(const char *name) const
Definition: TF1.h:373
TF1Parameters(const TF1Parameters &rhs)
Definition: TF1.h:68