Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
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 <string>
28#include <vector>
29#include "TFormula.h"
30#include "TMethodCall.h"
31#include "TAttLine.h"
32#include "TAttFill.h"
33#include "TAttMarker.h"
34#include "TF1AbsComposition.h"
35#include "TMath.h"
36#include "Math/Types.h"
37#include "Math/ParamFunctor.h"
38
39class TF1;
40class TH1;
41class TAxis;
42class TRandom;
43
44namespace ROOT {
45 namespace Fit {
46 class FitResult;
47 }
48}
49
51public:
52 TF1Parameters() {} // needed for the I/O
54 fParameters(std::vector<Double_t>(npar)),
55 fParNames(std::vector<std::string>(npar))
56 {
57 for (int i = 0; i < npar; ++i) {
58 fParNames[i] = std::string(TString::Format("p%d", i).Data());
59 }
60 }
61 // copy constructor
65 {}
66 // assignment
68 {
69 if (&rhs == this) return *this;
71 fParNames = rhs.fParNames;
72 return *this;
73 }
74 virtual ~TF1Parameters() {}
75
76 // getter methods
78 {
79 return (CheckIndex(iparam)) ? fParameters[iparam] : 0;
80 }
81 Double_t GetParameter(const char *name) const
82 {
84 }
85 const Double_t *GetParameters() const
86 {
87 return fParameters.data();
88 }
89 const std::vector<double> &ParamsVec() const
90 {
91 return fParameters;
92 }
93
94 Int_t GetParNumber(const char *name) const;
95
96 const char *GetParName(Int_t iparam) const
97 {
98 return (CheckIndex(iparam)) ? fParNames[iparam].c_str() : "";
99 }
100
101
102 // setter methods
104 {
105 if (!CheckIndex(iparam)) return;
106 fParameters[iparam] = value;
107 }
108 void SetParameters(const Double_t *params)
109 {
110 std::copy(params, params + fParameters.size(), fParameters.begin());
111 }
112 void SetParameters(Double_t p0, Double_t p1, Double_t p2 = 0, Double_t p3 = 0, Double_t p4 = 0,
113 Double_t p5 = 0, Double_t p6 = 0, Double_t p7 = 0, Double_t p8 = 0,
114 Double_t p9 = 0, Double_t p10 = 0);
115
116 void SetParameter(const char *name, Double_t value)
117 {
119 }
120 void SetParName(Int_t iparam, const char *name)
121 {
122 if (!CheckIndex(iparam)) return;
123 fParNames[iparam] = std::string(name);
124 }
125 void SetParNames(const char *name0 = "p0", const char *name1 = "p1", const char *name2 = "p2",
126 const char *name3 = "p3", const char *name4 = "p4", const char *name5 = "p5",
127 const char *name6 = "p6", const char *name7 = "p7", const char *name8 = "p8",
128 const char *name9 = "p9", const char *name10 = "p10");
129
130
131
132 ClassDef(TF1Parameters, 1) // The Parameters of a parameteric function
133private:
134
135 bool CheckIndex(Int_t i) const
136 {
137 return (i >= 0 && i < int(fParameters.size()));
138 }
139
140 std::vector<Double_t> fParameters; // parameter values
141 std::vector<std::string> fParNames; // parameter names
142};
143
144namespace ROOT {
145 namespace Internal {
146 /// %Internal class used by TF1 for defining
147 /// template specialization for different TF1 constructors
148 template<class Func>
149 struct TF1Builder {
150 static void Build(TF1 *f, Func func);
151 };
152
153 template<class Func>
154 struct TF1Builder<Func *> {
155 static void Build(TF1 *f, Func *func);
156 };
157
158 /// %Internal class used by TF1 for obtaining the type from a functor
159 /// out of the set of valid operator() signatures.
160 template<typename T>
162 };
163
164 template<typename F, typename T>
165 struct GetFunctorType<T(F::*)(const T *, const double *)> {
166 using type = T;
167 };
168
169 template<typename F, typename T>
170 struct GetFunctorType<T(F::*)(const T *, const double *) const> {
171 using type = T;
172 };
173
174 template<typename F, typename T>
175 struct GetFunctorType<T(F::*)(T *, double *)> {
176 using type = T;
177 };
178
179 template<typename F, typename T>
180 struct GetFunctorType<T(F::*)(T *, double *) const> {
181 using type = T;
182 };
183
184 /// %Internal class used by TF1 to get the right operator() signature
185 /// from a Functor with several ones.
186 template<typename T, typename F>
187 auto GetTheRightOp(T(F::*opPtr)(const T *, const double *)) -> decltype(opPtr)
188 {
189 return opPtr;
190 }
191
192 template<typename T, typename F>
193 auto GetTheRightOp(T(F::*opPtr)(const T *, const double *) const) -> decltype(opPtr)
194 {
195 return opPtr;
196 }
197
198 template<typename T, typename F>
199 auto GetTheRightOp(T(F::*opPtr)(T *, double *)) -> decltype(opPtr)
200 {
201 return opPtr;
202 }
203
204 template<typename T, typename F>
205 auto GetTheRightOp(T(F::*opPtr)(T *, double *) const) -> decltype(opPtr)
206 {
207 return opPtr;
208 }
209 }
210}
211
212
213class TF1 : public TNamed, public TAttLine, public TAttFill, public TAttMarker {
214
215 template<class Func>
217
218public:
219 /// Add to list behavior
220 enum class EAddToList {
221 kDefault,
222 kAdd,
223 kNo
224 };
225
226
229 virtual TF1FunctorPointer * Clone() const = 0;
230 };
231
232protected:
233
234 enum EFType {
235 kFormula = 0, ///< Formula functions which can be stored,
236 kPtrScalarFreeFcn, ///< Pointer to scalar free function,
237 kInterpreted, ///< Interpreted functions constructed by name,
238 kTemplVec, ///< Vectorized free functions or TemplScalar functors evaluating on vectorized parameters,
239 kTemplScalar, ///< TemplScalar functors evaluating on scalar parameters
241 }; // formula based on composition class (e.g. NSUM, CONV)
242
243 Double_t fXmin{-1111}; ///< Lower bounds for the range
244 Double_t fXmax{-1111}; ///< Upper bounds for the range
245 Int_t fNpar{}; ///< Number of parameters
246 Int_t fNdim{}; ///< Function dimension
247 Int_t fNpx{100}; ///< Number of points used for the graphical representation
249 Int_t fNpfits{}; ///< Number of points used in the fit
250 Int_t fNDF{}; ///< Number of degrees of freedom in the fit
251 Double_t fChisquare{}; ///< Function fit chisquare
252 Double_t fMinimum{-1111}; ///< Minimum value for plotting
253 Double_t fMaximum{-1111}; ///< Maximum value for plotting
254 std::vector<Double_t> fParErrors; ///< Array of errors of the fNpar parameters
255 std::vector<Double_t> fParMin; ///< Array of lower limits of the fNpar parameters
256 std::vector<Double_t> fParMax; ///< Array of upper limits of the fNpar parameters
257 std::vector<Double_t> fSave; ///< Array of fNsave function values
258 std::vector<Double_t> fIntegral; ///<! Integral of function binned on fNpx bins
259 std::vector<Double_t> fAlpha; ///<! Array alpha. for each bin in x the deconvolution r of fIntegral
260 std::vector<Double_t> fBeta; ///<! Array beta. is approximated by x = alpha +beta*r *gamma*r**2
261 std::vector<Double_t> fGamma; ///<! Array gamma.
262 TObject *fParent{nullptr}; ///<! Parent object hooking this function (if one)
263 TH1 *fHistogram{nullptr}; ///<! Pointer to histogram used for visualisation
264 std::unique_ptr<TMethodCall> fMethodCall; ///<! Pointer to MethodCall in case of interpreted function
265 Bool_t fNormalized{false}; ///< Normalization option (false by default)
266 Double_t fNormIntegral{}; ///< Integral of the function before being normalized
267 std::unique_ptr<TF1FunctorPointer> fFunctor; ///<! Functor object to wrap any C++ callable object
268 std::unique_ptr<TFormula> fFormula; ///< Pointer to TFormula in case when user define formula
269 std::unique_ptr<TF1Parameters> fParams; ///< Pointer to Function parameters object (exists only for not-formula functions)
270 std::unique_ptr<TF1AbsComposition> fComposition; ///< Pointer to composition (NSUM or CONV)
271
272 /// General constructor for TF1. Most of the other constructors delegate on it
273 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):
274 TNamed(name, name), TAttLine(), TAttFill(), TAttMarker(), fXmin(xmin), fXmax(xmax), fNpar(npar), fNdim(ndim),
275 fType(functionType), fParErrors(npar), fParMin(npar), fParMax(npar)
276 {
277 fParams.reset(params);
278 fFunctor.reset(functor);
279 DoInitialize(addToGlobList);
280 }
281
282private:
283 // NSUM parsing helper functions
284 void DefineNSUMTerm(TObjArray *newFuncs, TObjArray *coeffNames,
285 TString &fullFormula,
286 TString &formula, int termStart, int termEnd,
288 int TermCoeffLength(TString &term);
289
290protected:
291
292 template <class T>
295 TF1FunctorPointerImpl(const std::function<T(const T *f, const Double_t *param)> &func) : fImpl(func){};
297 TF1FunctorPointer * Clone() const override { return new TF1FunctorPointerImpl<T>(fImpl); }
299 };
300
301
302
303
304 static std::atomic<Bool_t> fgAbsValue; //use absolute value of function when computing integral
305 static Bool_t fgRejectPoint; //True if point must be rejected in a fit
306 static std::atomic<Bool_t> fgAddToGlobList; //True if we want to register the function in the global list
307 static TF1 *fgCurrent; //pointer to current function being processed
308
309
310 //void CreateFromFunctor(const char *name, Int_t npar, Int_t ndim = 1);
311 void DoInitialize(EAddToList addToGlobList);
312
314 // tabulate the cumulative function integral at fNpx points. Used by GetRandom
316
317 virtual Double_t GetMinMaxNDim(Double_t *x , Bool_t findmax, Double_t epsilon = 0, Int_t maxiter = 0) const;
318 virtual void GetRange(Double_t *xmin, Double_t *xmax) const;
320
321public:
322
323 // TF1 status bits
325 kNotGlobal = BIT(10), // don't register in global list of functions
326 kNotDraw = BIT(9) // don't draw the function when in a TH1
327 };
328
329 TF1();
330 TF1(const char *name, const char *formula, Double_t xmin = 0, Double_t xmax = 1, EAddToList addToGlobList = EAddToList::kDefault, bool vectorize = false);
331 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
332 TF1(const char *name, Double_t xmin, Double_t xmax, Int_t npar, Int_t ndim = 1, EAddToList addToGlobList = EAddToList::kDefault);
333 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);
334 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);
335
336 template <class T>
337 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):
338 TF1(EFType::kTemplScalar, name, xmin, xmax, npar, ndim, addToGlobList, new TF1Parameters(npar), new TF1FunctorPointerImpl<T>(fcn))
339 {
340 fType = std::is_same<T, double>::value ? TF1::EFType::kTemplScalar : TF1::EFType::kTemplVec;
341 }
342
343 ////////////////////////////////////////////////////////////////////////////////
344 /// Constructor using a pointer to function.
345 ///
346 /// \param[in] name object name
347 /// \param[in] fcn pointer to function
348 /// \param[in] xmin,xmax x axis limits
349 /// \param[in] npar is the number of free parameters used by the function
350 /// \param[in] ndim number of dimensions
351 /// \param[in] addToGlobList boolean marking if it should be added to global list
352 ///
353 /// This constructor creates a function of type C when invoked
354 /// with the normal C++ compiler.
355 ///
356 ///
357 /// \warning A function created with this constructor cannot be Cloned
358
359
360 template <class T>
361 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):
362 TF1(EFType::kTemplVec, name, xmin, xmax, npar, ndim, addToGlobList, new TF1Parameters(npar), new TF1FunctorPointerImpl<T>(fcn))
363 {}
364
365 // Constructors using functors (compiled mode only)
366 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);
367
368 // Template constructors from any C++ callable object, defining the operator() (double * , double *)
369 // and returning a double.
370 // The class name is not needed when using compile code, while it is required when using
371 // interpreted code via the specialized constructor with void *.
372 // An instance of the C++ function class or its pointer can both be used. The former is reccomended when using
373 // C++ compiled code, but if CINT compatibility is needed, then a pointer to the function class must be used.
374 // xmin and xmax specify the plotting range, npar is the number of parameters.
375 // See the tutorial math/exampleFunctor.C for an example of using this constructor
376 template <typename Func>
377 TF1(const char *name, Func f, Double_t xmin, Double_t xmax, Int_t npar, Int_t ndim = 1, EAddToList addToGlobList = EAddToList::kDefault) :
378 TF1(EFType::kTemplScalar, name, xmin, xmax, npar, ndim, addToGlobList)
379 {
380 //actual fType set in TF1Builder
382 }
383
384 // backward compatible interface
385 template <typename Func>
386 TF1(const char *name, Func f, Double_t xmin, Double_t xmax, Int_t npar, const char *, EAddToList addToGlobList = EAddToList::kDefault) :
387 TF1(EFType::kTemplScalar, name, xmin, xmax, npar, 1, addToGlobList, new TF1Parameters(npar))
388 {
390 }
391
392
393 // Template constructors from a pointer to any C++ class of type PtrObj with a specific member function of type
394 // MemFn.
395 // The member function must have the signature of (double * , double *) and returning a double.
396 // The class name and the method name are not needed when using compile code
397 // (the member function pointer is used in this case), while they are required when using interpreted
398 // code via the specialized constructor with void *.
399 // xmin and xmax specify the plotting range, npar is the number of parameters.
400 // See the tutorial math/exampleFunctor.C for an example of using this constructor
401 template <class PtrObj, typename MemFn>
402 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) :
403 TF1(EFType::kTemplScalar, name, xmin, xmax, npar, ndim, addToGlobList, new TF1Parameters(npar), new TF1FunctorPointerImpl<double>(ROOT::Math::ParamFunctor(p, memFn)))
404 {}
405
406 // backward compatible interface
407 template <class PtrObj, typename MemFn>
408 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) :
409 TF1(EFType::kTemplScalar, name, xmin, xmax, npar, 1, addToGlobList, new TF1Parameters(npar), new TF1FunctorPointerImpl<double>(ROOT::Math::ParamFunctor(p, memFn)))
410 {}
411
412 TF1(const TF1 &f1);
413 TF1 &operator=(const TF1 &rhs);
414 ~TF1() override;
415 virtual void AddParameter(const TString &name, Double_t value)
416 {
417 if (fFormula) fFormula->AddParameter(name, value);
418 }
419 // virtual void AddParameters(const pair<TString,Double_t> *pairs, Int_t size) { fFormula->AddParameters(pairs,size); }
420 // virtual void AddVariable(const TString &name, Double_t value = 0) { if (fFormula) fFormula->AddVariable(name,value); }
421 // virtual void AddVariables(const TString *vars, Int_t size) { if (fFormula) fFormula->AddVariables(vars,size); }
424 void Browse(TBrowser *b) override;
425 void Copy(TObject &f1) const override;
426 TObject *Clone(const char *newname = nullptr) const override;
427 virtual Double_t Derivative(Double_t x, Double_t *params = nullptr, Double_t epsilon = 0.001) const;
428 virtual Double_t Derivative2(Double_t x, Double_t *params = nullptr, Double_t epsilon = 0.001) const;
429 virtual Double_t Derivative3(Double_t x, Double_t *params = nullptr, Double_t epsilon = 0.001) const;
430 static Double_t DerivativeError();
431 Int_t DistancetoPrimitive(Int_t px, Int_t py) override;
432 void Draw(Option_t *option = "") override;
433 virtual TF1 *DrawCopy(Option_t *option = "") const;
434 virtual TObject *DrawDerivative(Option_t *option = "al"); // *MENU*
435 virtual TObject *DrawIntegral(Option_t *option = "al"); // *MENU*
436 virtual void DrawF1(Double_t xmin, Double_t xmax, Option_t *option = "");
437 virtual Double_t Eval(Double_t x, Double_t y = 0, Double_t z = 0, Double_t t = 0) const;
438 //template <class T> T Eval(T x, T y = 0, T z = 0, T t = 0) const;
439 virtual Double_t EvalPar(const Double_t *x, const Double_t *params = nullptr);
440 template <class T> T EvalPar(const T *x, const Double_t *params = nullptr);
441 virtual Double_t operator()(Double_t x, Double_t y = 0, Double_t z = 0, Double_t t = 0) const;
442 template <class T> T operator()(const T *x, const Double_t *params = nullptr);
443 void ExecuteEvent(Int_t event, Int_t px, Int_t py) override;
444 virtual void FixParameter(Int_t ipar, Double_t value);
446 {
447 return (fType == EFType::kTemplVec) || (fType == EFType::kFormula && fFormula && fFormula->IsVectorized());
448 }
450 {
451 return fChisquare;
452 }
453 virtual TH1 *GetHistogram() const;
455 {
457 }
459 {
460 return fFormula.get();
461 }
462 virtual const TFormula *GetFormula() const
463 {
464 return fFormula.get();
465 }
467 {
468 return (fFormula) ? fFormula->GetExpFormula(option) : TString();
469 }
470 virtual const TObject *GetLinearPart(Int_t i) const
471 {
472 return (fFormula) ? fFormula->GetLinearPart(i) : nullptr;
473 }
474 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;
475 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;
476 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;
477 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;
479 {
480 return fMaximum;
481 }
483 {
484 return fMinimum;
485 }
486 virtual Int_t GetNpar() const
487 {
488 return fNpar;
489 }
490 virtual Int_t GetNdim() const
491 {
492 return fNdim;
493 }
494 virtual Int_t GetNDF() const;
495 virtual Int_t GetNpx() const
496 {
497 return fNpx;
498 }
500 {
501 return fMethodCall.get();
502 }
503 virtual Int_t GetNumber() const
504 {
505 return (fFormula) ? fFormula->GetNumber() : 0;
506 }
507 virtual Int_t GetNumberFreeParameters() const;
508 virtual Int_t GetNumberFitPoints() const
509 {
510 return fNpfits;
511 }
512 char *GetObjectInfo(Int_t px, Int_t py) const override;
514 {
515 return fParent;
516 }
517 virtual Double_t GetParameter(Int_t ipar) const
518 {
519 return (fFormula) ? fFormula->GetParameter(ipar) : fParams->GetParameter(ipar);
520 }
521 virtual Double_t GetParameter(const TString &name) const
522 {
523 return (fFormula) ? fFormula->GetParameter(name) : fParams->GetParameter(name);
524 }
525 virtual Double_t *GetParameters() const
526 {
527 return (fFormula) ? fFormula->GetParameters() : const_cast<Double_t *>(fParams->GetParameters());
528 }
529 virtual void GetParameters(Double_t *params)
530 {
531 if (fFormula) fFormula->GetParameters(params);
532 else std::copy(fParams->ParamsVec().begin(), fParams->ParamsVec().end(), params);
533 }
534 virtual const char *GetParName(Int_t ipar) const
535 {
536 return (fFormula) ? fFormula->GetParName(ipar) : fParams->GetParName(ipar);
537 }
538 virtual Int_t GetParNumber(const char *name) const
539 {
540 return (fFormula) ? fFormula->GetParNumber(name) : fParams->GetParNumber(name);
541 }
542 virtual Double_t GetParError(Int_t ipar) const;
543 virtual const Double_t *GetParErrors() const
544 {
545 return fParErrors.data();
546 }
547 virtual void GetParLimits(Int_t ipar, Double_t &parmin, Double_t &parmax) const;
548 virtual Double_t GetProb() const;
549 virtual Int_t GetQuantiles(Int_t nprobSum, Double_t *q, const Double_t *probSum);
550 virtual Double_t GetRandom(TRandom * rng = nullptr, Option_t * opt = nullptr);
551 virtual Double_t GetRandom(Double_t xmin, Double_t xmax, TRandom * rng = nullptr, Option_t * opt = nullptr);
552 virtual void GetRange(Double_t &xmin, Double_t &xmax) const;
553 virtual void GetRange(Double_t &xmin, Double_t &ymin, Double_t &xmax, Double_t &ymax) const;
554 virtual void GetRange(Double_t &xmin, Double_t &ymin, Double_t &zmin, Double_t &xmax, Double_t &ymax, Double_t &zmax) const;
555 virtual Double_t GetSave(const Double_t *x);
556 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;
557 virtual Double_t GetXmin() const
558 {
559 return fXmin;
560 }
561 virtual Double_t GetXmax() const
562 {
563 return fXmax;
564 }
565 TAxis *GetXaxis() const ;
566 TAxis *GetYaxis() const ;
567 TAxis *GetZaxis() const ;
569 {
570 return (fFormula) ? fFormula->GetVariable(name) : 0;
571 }
572 virtual Double_t GradientPar(Int_t ipar, const Double_t *x, Double_t eps = 0.01);
573 template <class T>
574 T GradientPar(Int_t ipar, const T *x, Double_t eps = 0.01);
575 template <class T>
576 T GradientParTempl(Int_t ipar, const T *x, Double_t eps = 0.01);
577
578 virtual void GradientPar(const Double_t *x, Double_t *grad, Double_t eps = 0.01);
579 template <class T>
580 void GradientPar(const T *x, T *grad, Double_t eps = 0.01);
581 template <class T>
582 void GradientParTempl(const T *x, T *grad, Double_t eps = 0.01);
583
584 virtual void InitArgs(const Double_t *x, const Double_t *params);
585 static void InitStandardFunctions();
586 virtual Double_t Integral(Double_t a, Double_t b, Double_t epsrel = 1.e-12);
587 virtual Double_t IntegralOneDim(Double_t a, Double_t b, Double_t epsrel, Double_t epsabs, Double_t &err);
588 virtual Double_t IntegralError(Double_t a, Double_t b, const Double_t *params = nullptr, const Double_t *covmat = nullptr, Double_t epsilon = 1.E-2);
589 virtual Double_t IntegralError(Int_t n, const Double_t *a, const Double_t *b, const Double_t *params = nullptr, const Double_t *covmat = nullptr, Double_t epsilon = 1.E-2);
590 // virtual Double_t IntegralFast(const TGraph *g, Double_t a, Double_t b, Double_t *params = nullptr);
591 virtual Double_t IntegralFast(Int_t num, Double_t *x, Double_t *w, Double_t a, Double_t b, Double_t *params = nullptr, Double_t epsilon = 1e-12);
592 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);
593 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)
594 {
595 return IntegralMultiple(n, a, b, maxpts, epsrel, epsrel, relerr, nfnevl, ifail);
596 }
597 virtual Double_t IntegralMultiple(Int_t n, const Double_t *a, const Double_t *b, Double_t epsrel, Double_t &relerr);
598 virtual Bool_t IsEvalNormalized() const
599 {
600 return fNormalized;
601 }
602 /// return kTRUE if the point is inside the function range
603 virtual Bool_t IsInside(const Double_t *x) const
604 {
605 return !((x[0] < fXmin) || (x[0] > fXmax));
606 }
607 virtual Bool_t IsLinear() const
608 {
609 return (fFormula) ? fFormula->IsLinear() : false;
610 }
611 virtual Bool_t IsValid() const;
612 void Print(Option_t *option = "") const override;
613 void Paint(Option_t *option = "") override;
614 virtual void ReleaseParameter(Int_t ipar);
616 void SavePrimitive(std::ostream &out, Option_t *option = "") override;
617 virtual void SetChisquare(Double_t chi2)
618 {
619 fChisquare = chi2;
620 }
621 virtual void SetFitResult(const ROOT::Fit::FitResult &result, const Int_t *indpar = nullptr);
622 template <class PtrObj, typename MemFn>
623 void SetFunction(PtrObj &p, MemFn memFn);
624 template <typename Func>
625 void SetFunction(Func f);
626 virtual void SetMaximum(Double_t maximum = -1111); // *MENU*
627 virtual void SetMinimum(Double_t minimum = -1111); // *MENU*
628 virtual void SetNDF(Int_t ndf);
629 virtual void SetNumberFitPoints(Int_t npfits)
630 {
631 fNpfits = npfits;
632 }
633 virtual void SetNormalized(Bool_t flag)
634 {
635 fNormalized = flag;
636 Update();
637 }
638 virtual void SetNpx(Int_t npx = 100); // *MENU*
639 virtual void SetParameter(Int_t param, Double_t value)
640 {
641 (fFormula) ? fFormula->SetParameter(param, value) : fParams->SetParameter(param, value);
642 Update();
643 }
644 virtual void SetParameter(const TString &name, Double_t value)
645 {
646 (fFormula) ? fFormula->SetParameter(name, value) : fParams->SetParameter(name, value);
647 Update();
648 }
649 virtual void SetParameters(const Double_t *params)
650 {
651 (fFormula) ? fFormula->SetParameters(params) : fParams->SetParameters(params);
652 Update();
653 }
654 virtual void SetParameters(double p0, double p1 = 0.0, double p2 = 0.0, double p3 = 0.0, double p4 = 0.0,
655 double p5 = 0.0, double p6 = 0.0, double p7 = 0.0, double p8 = 0.0,
656 double p9 = 0.0, double p10 = 0.0)
657 {
658 if (fFormula) fFormula->SetParameters(p0, p1, p2, p3, p4, p5, p6, p7, p8, p9, p10);
659 else fParams->SetParameters(p0, p1, p2, p3, p4, p5, p6, p7, p8, p9, p10);
660 Update();
661 } // *MENU*
662 virtual void SetParName(Int_t ipar, const char *name);
663 virtual void SetParNames(const char *name0 = "p0", const char *name1 = "p1", const char *name2 = "p2",
664 const char *name3 = "p3", const char *name4 = "p4", const char *name5 = "p5",
665 const char *name6 = "p6", const char *name7 = "p7", const char *name8 = "p8",
666 const char *name9 = "p9", const char *name10 = "p10"); // *MENU*
667 virtual void SetParError(Int_t ipar, Double_t error);
668 virtual void SetParErrors(const Double_t *errors);
669 virtual void SetParLimits(Int_t ipar, Double_t parmin, Double_t parmax);
670 virtual void SetParent(TObject *p = nullptr)
671 {
672 fParent = p;
673 }
674 virtual void SetRange(Double_t xmin, Double_t xmax); // *MENU*
677 virtual void SetSavedPoint(Int_t point, Double_t value);
678 void SetTitle(const char *title = "") override; // *MENU*
679 virtual void SetVectorized(Bool_t vectorized)
680 {
682 fFormula->SetVectorized(vectorized);
683 else
684 Warning("SetVectorized", "Can only set vectorized flag on formula-based TF1");
685 }
686 virtual void Update();
687
688 static TF1 *GetCurrent();
689 static void AbsValue(Bool_t reject = kTRUE);
690 static void RejectPoint(Bool_t reject = kTRUE);
691 static Bool_t RejectedPoint();
692 static void SetCurrent(TF1 *f1);
693
694 //Moments
695 virtual Double_t Moment(Double_t n, Double_t a, Double_t b, const Double_t *params = nullptr, Double_t epsilon = 0.000001);
696 virtual Double_t CentralMoment(Double_t n, Double_t a, Double_t b, const Double_t *params = nullptr, Double_t epsilon = 0.000001);
697 virtual Double_t Mean(Double_t a, Double_t b, const Double_t *params = nullptr, Double_t epsilon = 0.000001)
698 {
699 return Moment(1, a, b, params, epsilon);
700 }
701 virtual Double_t Variance(Double_t a, Double_t b, const Double_t *params = nullptr, Double_t epsilon = 0.000001)
702 {
703 return CentralMoment(2, a, b, params, epsilon);
704 }
705
706 //some useful static utility functions to compute sampling points for Integral
707 //static void CalcGaussLegendreSamplingPoints(TGraph *g, Double_t eps=3.0e-11);
708 //static TGraph *CalcGaussLegendreSamplingPoints(Int_t num=21, Double_t eps=3.0e-11);
709 static void CalcGaussLegendreSamplingPoints(Int_t num, Double_t *x, Double_t *w, Double_t eps = 3.0e-11);
710
711private:
712 template <class T>
713 T EvalParTempl(const T *data, const Double_t *params = 0);
714
715#ifdef R__HAS_VECCORE
716 inline double EvalParVec(const Double_t *data, const Double_t *params);
717#endif
718
719 ClassDefOverride(TF1, 12) // The Parametric 1-D function
720};
721
722namespace ROOT {
723 namespace Internal {
724
725 template<class Func>
726 void TF1Builder<Func>::Build(TF1 *f, Func func)
727 {
728 using Fnc_t = typename ROOT::Internal::GetFunctorType<decltype(ROOT::Internal::GetTheRightOp(&Func::operator()))>::type;
729 f->fType = std::is_same<Fnc_t, double>::value? TF1::EFType::kTemplScalar : TF1::EFType::kTemplVec;
731 f->fParams.reset(new TF1Parameters(f->fNpar));
732 }
733
734 template<class Func>
735 void TF1Builder<Func *>::Build(TF1 *f, Func *func)
736 {
737 using Fnc_t = typename ROOT::Internal::GetFunctorType<decltype(ROOT::Internal::GetTheRightOp(&Func::operator()))>::type;
738 f->fType = std::is_same<Fnc_t, double>::value? TF1::EFType::kTemplScalar : TF1::EFType::kTemplVec;
740 f->fParams.reset(new TF1Parameters(f->fNpar));
741 }
742
743 /// TF1 building from a string
744 /// used to build a TFormula based on a lambda function
745 template<>
746 struct TF1Builder<const char *> {
747 static void Build(TF1 *f, const char *formula)
748 {
749 f->fType = TF1::EFType::kFormula;
750 f->fFormula.reset(new TFormula("tf1lambda", formula, f->fNdim, f->fNpar, false));
751 TString formulaExpression(formula);
752 Ssiz_t first = formulaExpression.Index("return") + 7;
753 Ssiz_t last = formulaExpression.Last(';');
754 TString title = formulaExpression(first, last - first);
755 f->SetTitle(title);
756 }
757 };
758 }
759}
760
762{
763 return Eval(x, y, z, t);
764}
765
766template <class T>
767inline T TF1::operator()(const T *x, const Double_t *params)
768{
769 return EvalPar(x, params);
770}
771
772////////////////////////////////////////////////////////////////////////////////
773/// EvalPar for vectorized
774template <class T>
775T TF1::EvalPar(const T *x, const Double_t *params)
776{
778 return EvalParTempl(x, params);
779 } else if (fType == EFType::kFormula) {
780 return fFormula->EvalPar(x, params);
781 } else
782 return TF1::EvalPar((double *)x, params);
783}
784
785////////////////////////////////////////////////////////////////////////////////
786/// Eval for vectorized functions
787// template <class T>
788// T TF1::Eval(T x, T y, T z, T t) const
789// {
790// if (fType == EFType::kFormula)
791// return fFormula->Eval(x, y, z, t);
792
793// T xx[] = {x, y, z, t};
794// Double_t *pp = (Double_t *)fParams->GetParameters();
795// return ((TF1 *)this)->EvalPar(xx, pp);
796// }
797
798// Internal to TF1. Evaluates Templated interfaces
799template <class T>
800inline T TF1::EvalParTempl(const T *data, const Double_t *params)
801{
803 if (!params) params = (Double_t *)fParams->GetParameters();
804 if (fFunctor)
805 return ((TF1FunctorPointerImpl<T> *)fFunctor.get())->fImpl(data, params);
806
807 // this should throw an error
808 // we nned to implement a vectorized GetSave(x)
809 return TMath::SignalingNaN();
810}
811
812#ifdef R__HAS_VECCORE
813// Internal to TF1. Evaluates Vectorized TF1 on data of type Double_v
814inline double TF1::EvalParVec(const Double_t *data, const Double_t *params)
815{
816 assert(fType == EFType::kTemplVec);
817 std::vector<ROOT::Double_v> d(fNdim);
818 ROOT::Double_v res;
819
820 for(auto i=0; i<fNdim; i++) {
821 d[i] = ROOT::Double_v(data[i]);
822 }
823
824 if (fFunctor) {
825 res = ((TF1FunctorPointerImpl<ROOT::Double_v> *) fFunctor.get())->fImpl(d.data(), params);
826 } else {
827 // res = GetSave(x);
828 return TMath::SignalingNaN();
829 }
830 return vecCore::Get<ROOT::Double_v>(res, 0);
831}
832#endif
833
835{
837}
839{
841}
842
843template <typename Func>
845{
846 // set function from a generic C++ callable object
849}
850template <class PtrObj, typename MemFn>
851void TF1::SetFunction(PtrObj &p, MemFn memFn)
852{
853 // set from a pointer to a member function
856}
857
858template <class T>
859inline T TF1::GradientPar(Int_t ipar, const T *x, Double_t eps)
860{
862 return GradientParTempl<T>(ipar, x, eps);
863 } else
864 return GradientParTempl<Double_t>(ipar, (const Double_t *)x, eps);
865}
866
867template <class T>
868inline T TF1::GradientParTempl(Int_t ipar, const T *x, Double_t eps)
869{
870 if (GetNpar() == 0)
871 return 0;
872
873 if (eps < 1e-10 || eps > 1) {
874 Warning("Derivative", "parameter esp=%g out of allowed range[1e-10,1], reset to 0.01", eps);
875 eps = 0.01;
876 }
877 Double_t h;
878 TF1 *func = (TF1 *)this;
879 Double_t *parameters = GetParameters();
880
881 // Copy parameters for thread safety
882 std::vector<Double_t> parametersCopy(parameters, parameters + GetNpar());
883 parameters = parametersCopy.data();
884
885 Double_t al, bl, h2;
886 T f1, f2, g1, g2, d0, d2;
887
888 ((TF1 *)this)->GetParLimits(ipar, al, bl);
889 if (al * bl != 0 && al >= bl) {
890 // this parameter is fixed
891 return 0;
892 }
893
894 // check if error has been computer (is not zero)
895 if (func->GetParError(ipar) != 0)
896 h = eps * func->GetParError(ipar);
897 else
898 h = eps;
899
900 // save original parameters
901 Double_t par0 = parameters[ipar];
902
903 parameters[ipar] = par0 + h;
904 f1 = func->EvalPar(x, parameters);
905 parameters[ipar] = par0 - h;
906 f2 = func->EvalPar(x, parameters);
907 parameters[ipar] = par0 + h / 2;
908 g1 = func->EvalPar(x, parameters);
909 parameters[ipar] = par0 - h / 2;
910 g2 = func->EvalPar(x, parameters);
911
912 // compute the central differences
913 h2 = 1 / (2. * h);
914 d0 = f1 - f2;
915 d2 = 2 * (g1 - g2);
916
917 T grad = h2 * (4 * d2 - d0) / 3.;
918
919 // restore original value
920 parameters[ipar] = par0;
921
922 return grad;
923}
924
925template <class T>
926inline void TF1::GradientPar(const T *x, T *grad, Double_t eps)
927{
929 GradientParTempl<T>(x, grad, eps);
930 } else
931 GradientParTempl<Double_t>((const Double_t *)x, (Double_t *)grad, eps);
932}
933
934template <class T>
935inline void TF1::GradientParTempl(const T *x, T *grad, Double_t eps)
936{
937 if (eps < 1e-10 || eps > 1) {
938 Warning("Derivative", "parameter esp=%g out of allowed range[1e-10,1], reset to 0.01", eps);
939 eps = 0.01;
940 }
941
942 for (Int_t ipar = 0; ipar < GetNpar(); ipar++) {
943 grad[ipar] = GradientParTempl<T>(ipar, x, eps);
944 }
945}
946
947#endif
#define d(i)
Definition RSha256.hxx:102
#define b(i)
Definition RSha256.hxx:100
#define f(i)
Definition RSha256.hxx:104
#define a(i)
Definition RSha256.hxx:99
#define h(i)
Definition RSha256.hxx:106
#define e(i)
Definition RSha256.hxx:103
bool Bool_t
Definition RtypesCore.h:63
int Int_t
Definition RtypesCore.h:45
constexpr Bool_t kFALSE
Definition RtypesCore.h:101
double Double_t
Definition RtypesCore.h:59
constexpr Bool_t kTRUE
Definition RtypesCore.h:100
const char Option_t
Definition RtypesCore.h:66
#define ClassDef(name, id)
Definition Rtypes.h:337
#define BIT(n)
Definition Rtypes.h:85
#define ClassDefOverride(name, id)
Definition Rtypes.h:341
winID h TVirtualViewer3D TVirtualGLPainter p
Option_t Option_t option
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void data
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void char Point_t Rectangle_t WindowAttributes_t Float_t Float_t Float_t Int_t Int_t UInt_t UInt_t Rectangle_t result
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void on
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void value
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void char Point_t Rectangle_t WindowAttributes_t Float_t Float_t Float_t Int_t Int_t UInt_t UInt_t Rectangle_t Int_t Int_t Window_t TString Int_t GCValues_t GetPrimarySelectionOwner GetDisplay GetScreen GetColormap GetNativeEvent const char const char dpyName wid window const char font_name cursor keysym reg const char only_if_exist regb h Point_t winding char text const char depth char const char Int_t count const char ColorStruct_t color const char Pixmap_t Pixmap_t PictureAttributes_t attr const char char ret_data h unsigned char height h Atom_t Int_t ULong_t ULong_t unsigned char prop_list Atom_t Atom_t Atom_t Time_t type
char name[80]
Definition TGX11.cxx:110
float xmin
float * q
float ymin
float xmax
float ymax
class containing the result of the fit and all the related information (fitted parameter values,...
Definition FitResult.h:47
Param Functor class for Multidimensional functions.
Fill Area Attributes class.
Definition TAttFill.h:19
Line Attributes class.
Definition TAttLine.h:18
Marker Attributes class.
Definition TAttMarker.h:19
Class to manage histogram axis.
Definition TAxis.h:30
Using a TBrowser one can browse all ROOT objects.
Definition TBrowser.h:37
TF1 Parameters class.
Definition TF1.h:50
std::vector< Double_t > fParameters
Definition TF1.h:140
const char * GetParName(Int_t iparam) const
Definition TF1.h:96
const std::vector< double > & ParamsVec() const
Definition TF1.h:89
virtual ~TF1Parameters()
Definition TF1.h:74
Double_t GetParameter(const char *name) const
Definition TF1.h:81
Double_t GetParameter(Int_t iparam) const
Definition TF1.h:77
void SetParName(Int_t iparam, const char *name)
Definition TF1.h:120
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:3841
TF1Parameters & operator=(const TF1Parameters &rhs)
Definition TF1.h:67
void SetParameter(Int_t iparam, Double_t value)
Definition TF1.h:103
void SetParameter(const char *name, Double_t value)
Definition TF1.h:116
TF1Parameters(Int_t npar)
Definition TF1.h:53
std::vector< std::string > fParNames
Definition TF1.h:141
TF1Parameters()
Definition TF1.h:52
bool CheckIndex(Int_t i) const
Definition TF1.h:135
void SetParameters(const Double_t *params)
Definition TF1.h:108
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:3809
TF1Parameters(const TF1Parameters &rhs)
Definition TF1.h:62
const Double_t * GetParameters() const
Definition TF1.h:85
1-Dim function class
Definition TF1.h:213
std::unique_ptr< TF1FunctorPointer > fFunctor
! Functor object to wrap any C++ callable object
Definition TF1.h:267
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,...
Definition TF1.cxx:1824
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:1697
virtual Double_t GetXmax() const
Definition TF1.h:561
virtual void ReleaseParameter(Int_t ipar)
Release parameter number ipar during a fit operation.
Definition TF1.cxx:3153
virtual void SetParError(Int_t ipar, Double_t error)
Set error for parameter number ipar.
Definition TF1.cxx:3471
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:3646
EAddToList
Add to list behavior.
Definition TF1.h:220
virtual Double_t Mean(Double_t a, Double_t b, const Double_t *params=nullptr, Double_t epsilon=0.000001)
Definition TF1.h:697
virtual Double_t Derivative(Double_t x, Double_t *params=nullptr, Double_t epsilon=0.001) const
Returns the first derivative of the function at point x, computed by Richardson's extrapolation metho...
Definition TF1.cxx:1114
virtual const Double_t * GetParErrors() const
Definition TF1.h:543
virtual Int_t GetNumber() const
Definition TF1.h:503
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:1890
std::vector< Double_t > fParErrors
Array of errors of the fNpar parameters.
Definition TF1.h:254
Int_t fNdim
Function dimension.
Definition TF1.h:246
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:3789
static void AbsValue(Bool_t reject=kTRUE)
Static function: set the fgAbsValue flag.
Definition TF1.cxx:985
virtual TH1 * GetHistogram() const
Return a pointer to the histogram used to visualise the function Note that this histogram is managed ...
Definition TF1.cxx:1585
virtual const TFormula * GetFormula() const
Definition TF1.h:462
virtual void GetParLimits(Int_t ipar, Double_t &parmin, Double_t &parmax) const
Return limits for parameter ipar.
Definition TF1.cxx:1941
Int_t fNpar
Number of parameters.
Definition TF1.h:245
virtual Double_t GetParameter(const TString &name) const
Definition TF1.h:521
TAxis * GetYaxis() const
Get y axis of the function.
Definition TF1.cxx:2410
virtual void GetParameters(Double_t *params)
Definition TF1.h:529
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:3417
virtual Double_t GetParError(Int_t ipar) const
Return value of parameter number ipar.
Definition TF1.cxx:1931
static std::atomic< Bool_t > fgAddToGlobList
Definition TF1.h:306
virtual Bool_t IsEvalNormalized() const
Definition TF1.h:598
virtual Double_t GetVariable(const TString &name)
Definition TF1.h:568
virtual Double_t IntegralError(Double_t a, Double_t b, const Double_t *params=nullptr, const Double_t *covmat=nullptr, Double_t epsilon=1.E-2)
Return Error on Integral of a parametric function between a and b due to the parameter uncertainties ...
Definition TF1.cxx:2706
virtual void SetChisquare(Double_t chi2)
Definition TF1.h:617
virtual Double_t IntegralFast(Int_t num, Double_t *x, Double_t *w, Double_t a, Double_t b, Double_t *params=nullptr, Double_t epsilon=1e-12)
Gauss-Legendre integral, see CalcGaussLegendreSamplingPoints.
Definition TF1.cxx:2776
Double_t fNormIntegral
Integral of the function before being normalized.
Definition TF1.h:266
Double_t GetChisquare() const
Definition TF1.h:449
virtual void SetMaximum(Double_t maximum=-1111)
Set the maximum value along Y for this function In case the function is already drawn,...
Definition TF1.cxx:3392
void Print(Option_t *option="") const override
This method must be overridden when a class wants to print itself.
Definition TF1.cxx:2895
T GradientParTempl(Int_t ipar, const T *x, Double_t eps=0.01)
Definition TF1.h:868
virtual TH1 * CreateHistogram()
Definition TF1.h:454
Double_t fXmin
Lower bounds for the range.
Definition TF1.h:243
std::unique_ptr< TMethodCall > fMethodCall
! Pointer to MethodCall in case of interpreted function
Definition TF1.h:264
virtual void Update()
Called by functions such as SetRange, SetNpx, SetParameters to force the deletion of the associated h...
Definition TF1.cxx:3611
virtual Double_t GetProb() const
Return the fit probability.
Definition TF1.cxx:1956
void IntegrateForNormalization()
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)
Definition TF1.h:337
virtual Double_t GradientPar(Int_t ipar, const Double_t *x, Double_t eps=0.01)
Compute the gradient (derivative) wrt a parameter ipar.
Definition TF1.cxx:2444
TAxis * GetZaxis() const
Get z axis of the function. (In case this object is a TF2 or TF3)
Definition TF1.cxx:2421
virtual Double_t GetRandom(TRandom *rng=nullptr, Option_t *opt=nullptr)
Return a random number following this function shape.
Definition TF1.cxx:2191
virtual void SetRange(Double_t xmin, Double_t xmax)
Initialize the upper and lower bounds to draw the function.
Definition TF1.cxx:3520
virtual Double_t GetMaximumStored() const
Definition TF1.h:478
virtual Int_t GetNpar() const
Definition TF1.h:486
virtual TString GetExpFormula(Option_t *option="") const
Definition TF1.h:466
std::vector< Double_t > fBeta
! Array beta. is approximated by x = alpha +beta*r *gamma*r**2
Definition TF1.h:260
Int_t fNDF
Number of degrees of freedom in the fit.
Definition TF1.h:250
TH1 * fHistogram
! Pointer to histogram used for visualisation
Definition TF1.h:263
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:593
std::unique_ptr< TF1AbsComposition > fComposition
Pointer to composition (NSUM or CONV)
Definition TF1.h:270
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:3482
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:3050
TObject * GetParent() const
Definition TF1.h:513
Int_t fNpfits
Number of points used in the fit.
Definition TF1.h:249
virtual Double_t Derivative2(Double_t x, Double_t *params=nullptr, Double_t epsilon=0.001) const
Returns the second derivative of the function at point x, computed by Richardson's extrapolation meth...
Definition TF1.cxx:1179
static void SetCurrent(TF1 *f1)
Static function setting the current function.
Definition TF1.cxx:3341
void SetFunction(PtrObj &p, MemFn memFn)
Definition TF1.h:851
virtual void SetParent(TObject *p=nullptr)
Definition TF1.h:670
std::vector< Double_t > fAlpha
! Array alpha. for each bin in x the deconvolution r of fIntegral
Definition TF1.h:259
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:408
virtual Double_t Integral(Double_t a, Double_t b, Double_t epsrel=1.e-12)
IntegralOneDim or analytical integral.
Definition TF1.cxx:2529
void SetTitle(const char *title="") override
Set function title if title has the form "fffffff;xxxx;yyyy", it is assumed that the function title i...
Definition TF1.cxx:3550
virtual Int_t GetNumberFitPoints() const
Definition TF1.h:508
std::unique_ptr< TFormula > fFormula
Pointer to TFormula in case when user define formula.
Definition TF1.h:268
static Double_t DerivativeError()
Static function returning the error of the last call to the of Derivative's functions.
Definition TF1.cxx:1278
std::vector< Double_t > fParMin
Array of lower limits of the fNpar parameters.
Definition TF1.h:255
static void InitStandardFunctions()
Create the basic function objects.
Definition TF1.cxx:2496
Double_t fMaximum
Maximum value for plotting.
Definition TF1.h:253
virtual void SetNpx(Int_t npx=100)
Set the number of points used to draw the function.
Definition TF1.cxx:3431
virtual Double_t * GetParameters() const
Definition TF1.h:525
Double_t fMinimum
Minimum value for plotting.
Definition TF1.h:252
int TermCoeffLength(TString &term)
Definition TF1.cxx:925
static Bool_t fgRejectPoint
Definition TF1.h:305
void Copy(TObject &f1) const override
Copy this F1 to a new F1.
Definition TF1.cxx:1006
virtual void SetNumberFitPoints(Int_t npfits)
Definition TF1.h:629
virtual Double_t GetMinimumStored() const
Definition TF1.h:482
virtual void SetNormalized(Bool_t flag)
Definition TF1.h:633
void Paint(Option_t *option="") override
Paint this function with its current attributes.
Definition TF1.cxx:2951
TF1 & operator=(const TF1 &rhs)
Operator =.
Definition TF1.cxx:943
virtual Int_t GetNumberFreeParameters() const
Return the number of free parameters.
Definition TF1.cxx:1901
virtual Double_t Moment(Double_t n, Double_t a, Double_t b, const Double_t *params=nullptr, Double_t epsilon=0.000001)
Return nth moment of function between a and b.
Definition TF1.cxx:3665
virtual Double_t CentralMoment(Double_t n, Double_t a, Double_t b, const Double_t *params=nullptr, 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:3702
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:402
Double_t fChisquare
Function fit chisquare.
Definition TF1.h:251
EStatusBits
Definition TF1.h:324
@ kNotGlobal
Definition TF1.h:325
@ kNotDraw
Definition TF1.h:326
virtual TFormula * GetFormula()
Definition TF1.h:458
virtual void InitArgs(const Double_t *x, const Double_t *params)
Initialize parameters addresses.
Definition TF1.cxx:2481
TF1(const char *name, Func f, Double_t xmin, Double_t xmax, Int_t npar, const char *, EAddToList addToGlobList=EAddToList::kDefault)
Definition TF1.h:386
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:2849
TMethodCall * GetMethodCall() const
Definition TF1.h:499
Int_t DistancetoPrimitive(Int_t px, Int_t py) override
Compute distance from point px,py to a function.
Definition TF1.cxx:1294
virtual Double_t operator()(Double_t x, Double_t y=0, Double_t z=0, Double_t t=0) const
Definition TF1.h:761
EFType fType
Definition TF1.h:248
Bool_t fNormalized
Normalization option (false by default)
Definition TF1.h:265
virtual Double_t Variance(Double_t a, Double_t b, const Double_t *params=nullptr, Double_t epsilon=0.000001)
Definition TF1.h:701
virtual void SetMinimum(Double_t minimum=-1111)
Set the minimum value along Y for this function In case the function is already drawn,...
Definition TF1.cxx:3405
virtual void AddParameter(const TString &name, Double_t value)
Definition TF1.h:415
virtual void GetRange(Double_t *xmin, Double_t *xmax) const
Return range of a generic N-D function.
Definition TF1.cxx:2280
void Browse(TBrowser *b) override
Browse.
Definition TF1.cxx:994
virtual const char * GetParName(Int_t ipar) const
Definition TF1.h:534
~TF1() override
TF1 default destructor.
Definition TF1.cxx:954
static TF1 * fgCurrent
Definition TF1.h:307
virtual Double_t EvalPar(const Double_t *x, const Double_t *params=nullptr)
Evaluate function with given coordinates and parameters.
Definition TF1.cxx:1469
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)
General constructor for TF1. Most of the other constructors delegate on it.
Definition TF1.h:273
Int_t fNpx
Number of points used for the graphical representation.
Definition TF1.h:247
virtual void SetParLimits(Int_t ipar, Double_t parmin, Double_t parmax)
Set lower and upper limits for parameter ipar.
Definition TF1.cxx:3499
void DoInitialize(EAddToList addToGlobList)
Common initialization of the TF1.
Definition TF1.cxx:803
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:1864
static TF1 * GetCurrent()
Static function returning the current function being processed.
Definition TF1.cxx:1570
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:1993
virtual void SetParName(Int_t ipar, const char *name)
Set name of parameter number ipar.
Definition TF1.cxx:3448
char * GetObjectInfo(Int_t px, Int_t py) const override
Redefines TObject::GetObjectInfo.
Definition TF1.cxx:1919
void ExecuteEvent(Int_t event, Int_t px, Int_t py) override
Execute action corresponding to one event.
Definition TF1.cxx:1537
virtual Double_t GetSave(const Double_t *x)
Get value corresponding to X in array of fSave values.
Definition TF1.cxx:2343
static std::atomic< Bool_t > fgAbsValue
Definition TF1.h:304
virtual Bool_t IsLinear() const
Definition TF1.h:607
TF1()
TF1 default constructor.
Definition TF1.cxx:488
virtual TF1 * DrawCopy(Option_t *option="") const
Draw a copy of this function with its current attributes.
Definition TF1.cxx:1364
std::vector< Double_t > fParMax
Array of upper limits of the fNpar parameters.
Definition TF1.h:256
bool IsVectorized()
Definition TF1.h:445
virtual void SetParameters(double p0, double p1=0.0, double p2=0.0, double p3=0.0, double p4=0.0, double p5=0.0, double p6=0.0, double p7=0.0, double p8=0.0, double p9=0.0, double p10=0.0)
Definition TF1.h:654
void SavePrimitive(std::ostream &out, Option_t *option="") override
Save primitive as a C++ statement(s) on output stream out.
Definition TF1.cxx:3216
virtual Bool_t IsValid() const
Return kTRUE if the function is valid.
Definition TF1.cxx:2880
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:839
std::vector< Double_t > fSave
Array of fNsave function values.
Definition TF1.h:257
static Bool_t RejectedPoint()
See TF1::RejectPoint above.
Definition TF1.cxx:3655
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:884
std::vector< Double_t > fGamma
! Array gamma.
Definition TF1.h:261
TObject * fParent
! Parent object hooking this function (if one)
Definition TF1.h:262
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:1724
virtual void DrawF1(Double_t xmin, Double_t xmax, Option_t *option="")
Draw function between xmin and xmax.
Definition TF1.cxx:1421
Bool_t ComputeCdfTable(Option_t *opt)
Compute the cumulative function at fNpx points between fXmin and fXmax.
Definition TF1.cxx:2080
virtual void SetParameters(const Double_t *params)
Definition TF1.h:649
virtual TObject * DrawIntegral(Option_t *option="al")
Draw integral of this function.
Definition TF1.cxx:1408
std::vector< Double_t > fIntegral
! Integral of function binned on fNpx bins
Definition TF1.h:258
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:3460
T EvalParTempl(const T *data, const Double_t *params=0)
Eval for vectorized functions.
Definition TF1.h:800
virtual TObject * DrawDerivative(Option_t *option="al")
Draw derivative of this function.
Definition TF1.cxx:1386
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:1440
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:1615
std::unique_ptr< TF1Parameters > fParams
Pointer to Function parameters object (exists only for not-formula functions)
Definition TF1.h:269
virtual void SetParameter(const TString &name, Double_t value)
Definition TF1.h:644
virtual void SetParameter(Int_t param, Double_t value)
Definition TF1.h:639
virtual Double_t Derivative3(Double_t x, Double_t *params=nullptr, Double_t epsilon=0.001) const
Returns the third derivative of the function at point x, computed by Richardson's extrapolation metho...
Definition TF1.cxx:1244
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:3163
TObject * Clone(const char *newname=nullptr) const override
Make a complete copy of the underlying object.
Definition TF1.cxx:1065
EFType
Definition TF1.h:234
@ kCompositionFcn
Definition TF1.h:240
@ kFormula
Formula functions which can be stored,.
Definition TF1.h:235
@ kPtrScalarFreeFcn
Pointer to scalar free function,.
Definition TF1.h:236
@ kTemplScalar
TemplScalar functors evaluating on scalar parameters.
Definition TF1.h:239
@ kTemplVec
Vectorized free functions or TemplScalar functors evaluating on vectorized parameters,...
Definition TF1.h:238
@ kInterpreted
Interpreted functions constructed by name,.
Definition TF1.h:237
virtual void SetSavedPoint(Int_t point, Double_t value)
Restore value of function saved at point.
Definition TF1.cxx:3534
virtual void FixParameter(Int_t ipar, Double_t value)
Fix the value of a parameter for a fit operation The specified value will be used in the fit and the ...
Definition TF1.cxx:1558
virtual Bool_t IsInside(const Double_t *x) const
return kTRUE if the point is inside the function range
Definition TF1.h:603
virtual Int_t GetNpx() const
Definition TF1.h:495
Double_t fXmax
Upper bounds for the range.
Definition TF1.h:244
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:1656
virtual Int_t GetNdim() const
Definition TF1.h:490
virtual Double_t GetXmin() const
Definition TF1.h:557
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:848
virtual const TObject * GetLinearPart(Int_t i) const
Definition TF1.h:470
virtual void SetVectorized(Bool_t vectorized)
Definition TF1.h:679
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:2619
TF1(const char *name, Func f, Double_t xmin, Double_t xmax, Int_t npar, Int_t ndim=1, EAddToList addToGlobList=EAddToList::kDefault)
Definition TF1.h:377
virtual Double_t GetParameter(Int_t ipar) const
Definition TF1.h:517
virtual Int_t GetParNumber(const char *name) const
Definition TF1.h:538
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)
Constructor using a pointer to function.
Definition TF1.h:361
virtual void SetFitResult(const ROOT::Fit::FitResult &result, const Int_t *indpar=nullptr)
Set the result from the fit parameter values, errors, chi2, etc... Optionally a pointer to a vector (...
Definition TF1.cxx:3353
TAxis * GetXaxis() const
Get x axis of the function.
Definition TF1.cxx:2399
The Formula class.
Definition TFormula.h:87
TH1 is the base class of all histogram classes in ROOT.
Definition TH1.h:58
Method or function calling interface.
Definition TMethodCall.h:37
The TNamed class is the base class for all named ROOT classes.
Definition TNamed.h:29
An array of TObjects.
Definition TObjArray.h:31
Mother of all ROOT objects.
Definition TObject.h:41
virtual void Warning(const char *method, const char *msgfmt,...) const
Issue warning message.
Definition TObject.cxx:956
This is the base class for the ROOT Random number generators.
Definition TRandom.h:27
Basic string class.
Definition TString.h:139
Ssiz_t Last(char c) const
Find last occurrence of a character c.
Definition TString.cxx:924
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:2356
Ssiz_t Index(const char *pat, Ssiz_t i=0, ECaseCompare cmp=kExact) const
Definition TString.h:651
Double_t y[n]
Definition legend1.C:17
Double_t x[n]
Definition legend1.C:17
const Int_t n
Definition legend1.C:16
TF1 * f1
Definition legend1.C:11
#define F(x, y, z)
TFitResultPtr Fit(FitObject *h1, TF1 *f1, Foption_t &option, const ROOT::Math::MinimizerOptions &moption, const char *goption, ROOT::Fit::DataRange &range)
Definition HFitImpl.cxx:133
Namespace for new Math classes and functions.
auto GetTheRightOp(T(F::*opPtr)(const T *, const double *)) -> decltype(opPtr)
Internal class used by TF1 to get the right operator() signature from a Functor with several ones.
Definition TF1.h:187
This file contains a specialised ROOT message handler to test for diagnostic in unit tests.
Double_t Double_v
Definition Types.h:56
Double_t SignalingNaN()
Returns a signaling NaN as defined by IEEE 754](http://en.wikipedia.org/wiki/NaN#Signaling_NaN).
Definition TMath.h:908
Definition first.py:1
Internal class used by TF1 for obtaining the type from a functor out of the set of valid operator() s...
Definition TF1.h:161
static void Build(TF1 *f, const char *formula)
Definition TF1.h:747
Internal class used by TF1 for defining template specialization for different TF1 constructors
Definition TF1.h:149
static void Build(TF1 *f, Func func)
Definition TF1.h:726
TF1FunctorPointer * Clone() const override
Definition TF1.h:297
ROOT::Math::ParamFunctorTempl< T > fImpl
Definition TF1.h:298
TF1FunctorPointerImpl(const std::function< T(const T *f, const Double_t *param)> &func)
Definition TF1.h:295
~TF1FunctorPointerImpl() override
Definition TF1.h:296
TF1FunctorPointerImpl(const ROOT::Math::ParamFunctorTempl< T > &func)
Definition TF1.h:294
virtual ~TF1FunctorPointer()
Definition TF1.h:228
virtual TF1FunctorPointer * Clone() const =0
th1 Draw()
double epsilon
Definition triangle.c:618