include/root/TF1.h

0001 // @(#)root/hist:$Id$
0002 // Author: Rene Brun   18/08/95
0003
0004 /*************************************************************************
0005  * Copyright (C) 1995-2000, Rene Brun and Fons Rademakers.               *
0006  * All rights reserved.                                                  *
0007  *                                                                       *
0008  * For the licensing terms see $ROOTSYS/LICENSE.                         *
0009  * For the list of contributors see $ROOTSYS/README/CREDITS.             *
0010  *************************************************************************/
0011 // ---------------------------------- F1.h
0012
0013 #ifndef ROOT_TF1
0014 #define ROOT_TF1
0015
0016 //////////////////////////////////////////////////////////////////////////
0017 //                                                                      //
0018 // TF1                                                                  //
0019 //                                                                      //
0020 // The Parametric 1-D function                                          //
0021 //                                                                      //
0022 //////////////////////////////////////////////////////////////////////////
0023
0024 #include "RConfigure.h"
0025 #include <functional>
0026 #include <cassert>
0027 #include <memory>
0028 #include <string>
0029 #include <vector>
0030 #include "TFormula.h"
0031 #include "TMethodCall.h"
0032 #include "TAttLine.h"
0033 #include "TAttFill.h"
0034 #include "TAttMarker.h"
0035 #include "TF1AbsComposition.h"
0036 #include "TMath.h"
0037 #include "Math/Types.h"
0038 #include "Math/ParamFunctor.h"
0039
0040 #include <string>
0041
0042 class TF1;
0043 class TH1;
0044 class TAxis;
0045 class TRandom;
0046
0047 namespace ROOT {
0048    namespace Fit {
0049       class FitResult;
0050    }
0051 }
0052
0053 class TF1Parameters {
0054 public:
0055    TF1Parameters() {} // needed for the I/O
0056    TF1Parameters(Int_t npar) :
0057       fParameters(std::vector<Double_t>(npar)),
0058       fParNames(std::vector<std::string>(npar))
0059    {
0060       for (int i = 0; i < npar; ++i) {
0061          fParNames[i] = std::string(TString::Format("p%d", i).Data());
0062       }
0063    }
0064    // copy constructor
0065    TF1Parameters(const TF1Parameters &rhs) :
0066       fParameters(rhs.fParameters),
0067       fParNames(rhs.fParNames)
0068    {}
0069    // assignment
0070    TF1Parameters &operator=(const TF1Parameters &rhs)
0071    {
0072       if (&rhs == this) return *this;
0073       fParameters = rhs.fParameters;
0074       fParNames = rhs.fParNames;
0075       return *this;
0076    }
0077    virtual ~TF1Parameters() {}
0078
0079    // getter methods
0080    Double_t GetParameter(Int_t iparam) const
0081    {
0082       return (CheckIndex(iparam)) ? fParameters[iparam] : 0;
0083    }
0084    Double_t GetParameter(const char *name) const
0085    {
0086       return GetParameter(GetParNumber(name));
0087    }
0088    const Double_t *GetParameters() const
0089    {
0090       return fParameters.data();
0091    }
0092    const std::vector<double> &ParamsVec() const
0093    {
0094       return fParameters;
0095    }
0096
0097    Int_t GetParNumber(const char *name) const;
0098
0099    const char *GetParName(Int_t iparam) const
0100    {
0101       return (CheckIndex(iparam)) ? fParNames[iparam].c_str() : "";
0102    }
0103
0104
0105    // setter methods
0106    void   SetParameter(Int_t iparam, Double_t value)
0107    {
0108       if (!CheckIndex(iparam)) return;
0109       fParameters[iparam] = value;
0110    }
0111    void  SetParameters(const Double_t *params)
0112    {
0113       std::copy(params, params + fParameters.size(), fParameters.begin());
0114    }
0115    template <typename... Args>
0116    void   SetParameters(Double_t arg1, Args &&... args);
0117
0118    void   SetParameter(const char *name, Double_t value)
0119    {
0120       SetParameter(GetParNumber(name), value);
0121    }
0122    void   SetParName(Int_t iparam, const char *name)
0123    {
0124       if (!CheckIndex(iparam)) return;
0125       fParNames[iparam] = std::string(name);
0126    }
0127    template <typename... Args>
0128    void   SetParNames(Args &&... args);
0129
0130    ClassDef(TF1Parameters, 1)  // The Parameters of a parameteric function
0131 private:
0132
0133    bool CheckIndex(Int_t i) const
0134    {
0135       return (i >= 0 && i < int(fParameters.size()));
0136    }
0137
0138    std::vector<Double_t> fParameters;    // parameter values
0139    std::vector<std::string> fParNames;   // parameter names
0140 };
0141
0142 /// Set parameter values.
0143 /// NaN values will be skipped, meaning that the corresponding parameters will not be changed.
0144 template <typename... Args>
0145 void TF1Parameters::SetParameters(Double_t arg1, Args &&...args)
0146 {
0147    int i = 0;
0148    for (double val : {arg1, static_cast<Double_t>(args)...}) {
0149       if(!TMath::IsNaN(val)) SetParameter(i++, val);
0150    }
0151 }
0152
0153 /// Set parameter names.
0154 /// Empty strings will be skipped, meaning that the corresponding name will not be changed.
0155 template <typename... Args>
0156 void TF1Parameters::SetParNames(Args &&...args)
0157 {
0158    int i = 0;
0159    for (auto name : {static_cast<std::string const&>(args)...}) {
0160       if(!name.empty()) SetParName(i++, name.c_str());
0161    }
0162 }
0163
0164 namespace ROOT {
0165    namespace Internal {
0166       /// %Internal class used by TF1 for defining
0167       /// template specialization for different TF1 constructors
0168       template<class Func>
0169       struct TF1Builder {
0170          static void Build(TF1 *f, Func func);
0171       };
0172
0173       template<class Func>
0174       struct TF1Builder<Func *> {
0175          static void Build(TF1 *f, Func *func);
0176       };
0177
0178       /// %Internal class used by TF1 for obtaining the type from a functor
0179       /// out of the set of valid operator() signatures.
0180       template<typename T>
0181       struct GetFunctorType {
0182       };
0183
0184       template<typename F, typename T>
0185       struct GetFunctorType<T(F::*)(const T *, const double *)> {
0186          using type = T;
0187       };
0188
0189       template<typename F, typename T>
0190       struct GetFunctorType<T(F::*)(const T *, const double *) const> {
0191          using type = T;
0192       };
0193
0194       template<typename F, typename T>
0195       struct GetFunctorType<T(F::*)(T *, double *)> {
0196          using type = T;
0197       };
0198
0199       template<typename F, typename T>
0200       struct GetFunctorType<T(F::*)(T *, double *) const> {
0201          using type = T;
0202       };
0203
0204       /// %Internal class used by TF1 to get the right operator() signature
0205       /// from a Functor with several ones.
0206       template<typename T, typename F>
0207       auto GetTheRightOp(T(F::*opPtr)(const T *, const double *)) -> decltype(opPtr)
0208       {
0209          return opPtr;
0210       }
0211
0212       template<typename T, typename F>
0213       auto GetTheRightOp(T(F::*opPtr)(const T *, const double *) const) -> decltype(opPtr)
0214       {
0215          return opPtr;
0216       }
0217
0218       template<typename T, typename F>
0219       auto GetTheRightOp(T(F::*opPtr)(T *, double *)) -> decltype(opPtr)
0220       {
0221          return opPtr;
0222       }
0223
0224       template<typename T, typename F>
0225       auto GetTheRightOp(T(F::*opPtr)(T *, double *) const) -> decltype(opPtr)
0226       {
0227          return opPtr;
0228       }
0229    }
0230 }
0231
0232
0233 class TF1 : public TNamed, public TAttLine, public TAttFill, public TAttMarker {
0234
0235    template<class Func>
0236    friend struct ROOT::Internal::TF1Builder;
0237
0238 public:
0239    /// Add to list behavior
0240    enum class EAddToList {
0241       kDefault,
0242       kAdd,
0243       kNo
0244    };
0245
0246
0247    struct TF1FunctorPointer {
0248       virtual  ~TF1FunctorPointer() {}
0249       virtual  TF1FunctorPointer * Clone() const = 0;
0250    };
0251
0252 protected:
0253
0254    enum EFType {
0255       kFormula = 0,      ///< Formula functions which can be stored,
0256       kPtrScalarFreeFcn, ///< Pointer to scalar free function,
0257       kInterpreted,      ///< Interpreted functions constructed by name,
0258       kTemplVec,         ///< Vectorized free functions or TemplScalar functors evaluating on vectorized parameters,
0259       kTemplScalar,      ///< TemplScalar functors evaluating on scalar parameters
0260       kCompositionFcn
0261    }; // formula based on composition class (e.g. NSUM, CONV)
0262
0263    Double_t    fXmin{-1111};                        ///<  Lower bounds for the range
0264    Double_t    fXmax{-1111};                        ///<  Upper bounds for the range
0265    Int_t       fNpar{};                             ///<  Number of parameters
0266    Int_t       fNdim{};                             ///<  Function dimension
0267    Int_t       fNpx{100};                           ///<  Number of points used for the graphical representation
0268    EFType      fType{EFType::kTemplScalar};
0269    Int_t       fNpfits{};                           ///<  Number of points used in the fit
0270    Int_t       fNDF{};                              ///<  Number of degrees of freedom in the fit
0271    Double_t    fChisquare{};                        ///<  Function fit chisquare
0272    Double_t    fMinimum{-1111};                     ///<  Minimum value for plotting
0273    Double_t    fMaximum{-1111};                     ///<  Maximum value for plotting
0274    std::vector<Double_t>    fParErrors;             ///<  Array of errors of the fNpar parameters
0275    std::vector<Double_t>    fParMin;                ///<  Array of lower limits of the fNpar parameters
0276    std::vector<Double_t>    fParMax;                ///<  Array of upper limits of the fNpar parameters
0277    std::vector<Double_t>    fSave;                  ///<  Array of fNsave function values
0278    std::vector<Double_t>    fIntegral;              ///<! Integral of function binned on fNpx bins
0279    std::vector<Double_t>    fAlpha;                 ///<! Array alpha. for each bin in x the deconvolution r of fIntegral
0280    std::vector<Double_t>    fBeta;                  ///<! Array beta.  is approximated by x = alpha +beta*r *gamma*r**2
0281    std::vector<Double_t>    fGamma;                 ///<! Array gamma.
0282    TObject     *fParent{nullptr};                   ///<! Parent object hooking this function (if one)
0283    TH1         *fHistogram{nullptr};                ///<! Pointer to histogram used for visualisation
0284    std::unique_ptr<TMethodCall> fMethodCall;        ///<! Pointer to MethodCall in case of interpreted function
0285    Bool_t      fNormalized{false};                  ///<  Normalization option (false by default)
0286    Double_t    fNormIntegral{};                     ///<  Integral of the function before being normalized
0287    std::unique_ptr<TF1FunctorPointer>  fFunctor;    ///<! Functor object to wrap any C++ callable object
0288    std::unique_ptr<TFormula>   fFormula;            ///<  Pointer to TFormula in case when user define formula
0289    std::unique_ptr<TF1Parameters> fParams;          ///<  Pointer to Function parameters object (exists only for not-formula functions)
0290    std::unique_ptr<TF1AbsComposition> fComposition; ///<  Pointer to composition (NSUM or CONV)
0291
0292    /// General constructor for TF1. Most of the other constructors delegate on it
0293    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):
0294       TNamed(name, name), TAttLine(), TAttFill(), TAttMarker(), fXmin(xmin), fXmax(xmax), fNpar(npar), fNdim(ndim),
0295       fType(functionType), fParErrors(npar), fParMin(npar), fParMax(npar)
0296    {
0297       fParams.reset(params);
0298       fFunctor.reset(functor);
0299       DoInitialize(addToGlobList);
0300    }
0301
0302 private:
0303    // NSUM parsing helper functions
0304    void DefineNSUMTerm(TObjArray *newFuncs, TObjArray *coeffNames,
0305                TString &fullFormula,
0306                TString &formula, int termStart, int termEnd,
0307                Double_t xmin, Double_t xmax);
0308    int TermCoeffLength(TString &term);
0309
0310 protected:
0311
0312    template <class T>
0313    struct TF1FunctorPointerImpl: TF1FunctorPointer {
0314       TF1FunctorPointerImpl(const ROOT::Math::ParamFunctorTempl<T> &func): fImpl(func) {};
0315       TF1FunctorPointerImpl(const std::function<T(const T *f, const Double_t *param)> &func) : fImpl(func){};
0316       ~TF1FunctorPointerImpl() override {}
0317        TF1FunctorPointer * Clone() const override { return new TF1FunctorPointerImpl<T>(fImpl); }
0318       ROOT::Math::ParamFunctorTempl<T> fImpl;
0319    };
0320
0321
0322
0323
0324    static std::atomic<Bool_t> fgAbsValue;  //use absolute value of function when computing integral
0325    static Bool_t fgRejectPoint;  //True if point must be rejected in a fit
0326    static std::atomic<Bool_t> fgAddToGlobList; //True if we want to register the function in the global list
0327    static TF1   *fgCurrent;   //pointer to current function being processed
0328
0329
0330    //void CreateFromFunctor(const char *name, Int_t npar, Int_t ndim = 1);
0331    void DoInitialize(EAddToList addToGlobList);
0332
0333    void IntegrateForNormalization();
0334    // tabulate the cumulative function integral at  fNpx points. Used by GetRandom
0335    Bool_t ComputeCdfTable(Option_t * opt);
0336
0337    virtual Double_t GetMinMaxNDim(Double_t *x , Bool_t findmax, Double_t epsilon = 0, Int_t maxiter = 0) const;
0338    virtual void GetRange(Double_t *xmin, Double_t *xmax) const;
0339    virtual TH1 *DoCreateHistogram(Double_t xmin, Double_t xmax, Bool_t recreate = kFALSE);
0340
0341 public:
0342
0343    // TF1 status bits
0344    enum EStatusBits {
0345       kNotGlobal   = BIT(10),  // don't register in global list of functions
0346       kNotDraw     = BIT(9)  // don't draw the function when in a TH1
0347    };
0348
0349    TF1();
0350    TF1(const char *name, const char *formula, Double_t xmin = 0, Double_t xmax = 1, EAddToList addToGlobList = EAddToList::kDefault, bool vectorize = false);
0351    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
0352    TF1(const char *name, Double_t xmin, Double_t xmax, Int_t npar, Int_t ndim = 1, EAddToList addToGlobList = EAddToList::kDefault);
0353    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);
0354    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);
0355
0356    template <class T>
0357    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):
0358       TF1(EFType::kTemplScalar, name, xmin, xmax, npar, ndim, addToGlobList, new TF1Parameters(npar), new TF1FunctorPointerImpl<T>(fcn))
0359    {
0360       fType = std::is_same<T, double>::value ? TF1::EFType::kTemplScalar : TF1::EFType::kTemplVec;
0361    }
0362
0363    ////////////////////////////////////////////////////////////////////////////////
0364    /// Constructor using a pointer to function.
0365    ///
0366    /// \param[in] name object name
0367    /// \param[in] fcn pointer to function
0368    /// \param[in] xmin,xmax x axis limits
0369    /// \param[in] npar is the number of free parameters used by the function
0370    /// \param[in] ndim number of dimensions
0371    /// \param[in] addToGlobList boolean marking if it should be added to global list
0372    ///
0373    /// This constructor creates a function of type C when invoked
0374    /// with the normal C++ compiler.
0375    ///
0376    ///
0377    /// \warning A function created with this constructor cannot be Cloned
0378
0379
0380    template <class T>
0381    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):
0382       TF1(EFType::kTemplVec, name, xmin, xmax, npar, ndim, addToGlobList, new TF1Parameters(npar), new TF1FunctorPointerImpl<T>(fcn))
0383    {}
0384
0385    // Constructors using functors (compiled mode only)
0386    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);
0387
0388    // Template constructors from any  C++ callable object,  defining  the operator() (double * , double *)
0389    // and returning a double.
0390    // The class name is not needed when using compile code, while it is required when using
0391    // interpreted code via the specialized constructor with void *.
0392    // An instance of the C++ function class or its pointer can both be used. The former is reccomended when using
0393    // C++ compiled code, but if CINT compatibility is needed, then a pointer to the function class must be used.
0394    // xmin and xmax specify the plotting range,  npar is the number of parameters.
0395    // See the tutorial math/exampleFunctor.C for an example of using this constructor
0396    template <typename Func>
0397    TF1(const char *name, Func f, Double_t xmin, Double_t xmax, Int_t npar, Int_t ndim = 1, EAddToList addToGlobList = EAddToList::kDefault) :
0398       TF1(EFType::kTemplScalar, name, xmin, xmax, npar, ndim, addToGlobList)
0399    {
0400       //actual fType set in TF1Builder
0401       ROOT::Internal::TF1Builder<Func>::Build(this, f);
0402    }
0403
0404    // backward compatible interface
0405    template <typename Func>
0406    TF1(const char *name, Func f, Double_t xmin, Double_t xmax, Int_t npar, const char *, EAddToList addToGlobList = EAddToList::kDefault) :
0407       TF1(EFType::kTemplScalar, name, xmin, xmax, npar, 1, addToGlobList, new TF1Parameters(npar))
0408    {
0409       ROOT::Internal::TF1Builder<Func>::Build(this, f);
0410    }
0411
0412
0413    // Template constructors from a pointer to any C++ class of type PtrObj with a specific member function of type
0414    // MemFn.
0415    // The member function must have the signature of  (double * , double *) and returning a double.
0416    // The class name and the method name are not needed when using compile code
0417    // (the member function pointer is used in this case), while they are required when using interpreted
0418    // code via the specialized constructor with void *.
0419    // xmin and xmax specify the plotting range,  npar is the number of parameters.
0420    // See the tutorial math/exampleFunctor.C for an example of using this constructor
0421    template <class PtrObj, typename MemFn>
0422    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) :
0423       TF1(EFType::kTemplScalar, name, xmin, xmax, npar, ndim, addToGlobList, new TF1Parameters(npar), new TF1FunctorPointerImpl<double>(ROOT::Math::ParamFunctor(p, memFn)))
0424    {}
0425
0426    // backward compatible interface
0427    template <class PtrObj, typename MemFn>
0428    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) :
0429       TF1(EFType::kTemplScalar, name, xmin, xmax, npar, 1, addToGlobList, new TF1Parameters(npar), new TF1FunctorPointerImpl<double>(ROOT::Math::ParamFunctor(p, memFn)))
0430    {}
0431
0432    TF1(const TF1 &f1);
0433    TF1 &operator=(const TF1 &rhs);
0434      ~TF1() override;
0435    virtual void     AddParameter(const TString &name, Double_t value)
0436    {
0437       if (fFormula) fFormula->AddParameter(name, value);
0438    }
0439    // virtual void     AddParameters(const pair<TString,Double_t> *pairs, Int_t size) { fFormula->AddParameters(pairs,size); }
0440    // virtual void     AddVariable(const TString &name, Double_t value = 0) { if (fFormula) fFormula->AddVariable(name,value); }
0441    // virtual void     AddVariables(const TString *vars, Int_t size) { if (fFormula) fFormula->AddVariables(vars,size); }
0442    virtual Bool_t   AddToGlobalList(Bool_t on = kTRUE);
0443    static  Bool_t   DefaultAddToGlobalList(Bool_t on = kTRUE);
0444    void     Browse(TBrowser *b) override;
0445    void     Copy(TObject &f1) const override;
0446    TObject         *Clone(const char *newname = nullptr) const override;
0447    virtual Double_t Derivative(Double_t x, Double_t *params = nullptr, Double_t epsilon = 0.001) const;
0448    virtual Double_t Derivative2(Double_t x, Double_t *params = nullptr, Double_t epsilon = 0.001) const;
0449    virtual Double_t Derivative3(Double_t x, Double_t *params = nullptr, Double_t epsilon = 0.001) const;
0450    static  Double_t DerivativeError();
0451    Int_t    DistancetoPrimitive(Int_t px, Int_t py) override;
0452    void     Draw(Option_t *option = "") override;
0453    virtual TF1     *DrawCopy(Option_t *option = "") const;
0454    virtual TObject *DrawDerivative(Option_t *option = "al"); // *MENU*
0455    virtual TObject *DrawIntegral(Option_t *option = "al"); // *MENU*
0456    virtual void     DrawF1(Double_t xmin, Double_t xmax, Option_t *option = "");
0457    virtual Double_t Eval(Double_t x, Double_t y = 0, Double_t z = 0, Double_t t = 0) const;
0458    //template <class T> T Eval(T x, T y = 0, T z = 0, T t = 0) const;
0459    virtual Double_t EvalPar(const Double_t *x, const Double_t *params = nullptr);
0460    template <class T> T EvalPar(const T *x, const Double_t *params = nullptr);
0461    virtual Double_t operator()(Double_t x, Double_t y = 0, Double_t z = 0, Double_t t = 0) const;
0462    template <class T> T operator()(const T *x, const Double_t *params = nullptr);
0463    void     ExecuteEvent(Int_t event, Int_t px, Int_t py) override;
0464    virtual void     FixParameter(Int_t ipar, Double_t value);
0465    bool      IsVectorized()
0466    {
0467       return (fType == EFType::kTemplVec) || (fType == EFType::kFormula && fFormula && fFormula->IsVectorized());
0468    }
0469    /// Return the Chisquare after fitting. See ROOT::Fit::FitResult::Chi2()
0470    Double_t     GetChisquare() const
0471    {
0472       return fChisquare;
0473    }
0474    virtual TH1     *GetHistogram() const;
0475    virtual TH1     *CreateHistogram()
0476    {
0477       return DoCreateHistogram(fXmin, fXmax);
0478    }
0479    virtual TFormula *GetFormula()
0480    {
0481       return fFormula.get();
0482    }
0483    virtual const TFormula *GetFormula() const
0484    {
0485       return fFormula.get();
0486    }
0487    virtual TString  GetExpFormula(Option_t *option = "") const
0488    {
0489       return (fFormula) ? fFormula->GetExpFormula(option) : TString();
0490    }
0491    virtual const TObject *GetLinearPart(Int_t i) const
0492    {
0493       return (fFormula) ? fFormula->GetLinearPart(i) : nullptr;
0494    }
0495    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;
0496    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;
0497    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;
0498    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;
0499    virtual Double_t GetMaximumStored() const
0500    {
0501       return fMaximum;
0502    }
0503    virtual Double_t GetMinimumStored() const
0504    {
0505       return fMinimum;
0506    }
0507    virtual Int_t    GetNpar() const
0508    {
0509       return fNpar;
0510    }
0511    virtual Int_t    GetNdim() const
0512    {
0513       return fNdim;
0514    }
0515    virtual Int_t    GetNDF() const;
0516    virtual Int_t    GetNpx() const
0517    {
0518       return fNpx;
0519    }
0520    TMethodCall    *GetMethodCall() const
0521    {
0522       return fMethodCall.get();
0523    }
0524    virtual Int_t    GetNumber() const
0525    {
0526       return (fFormula) ? fFormula->GetNumber() : 0;
0527    }
0528    virtual Int_t    GetNumberFreeParameters() const;
0529    virtual Int_t    GetNumberFitPoints() const
0530    {
0531       return fNpfits;
0532    }
0533    char    *GetObjectInfo(Int_t px, Int_t py) const override;
0534    TObject    *GetParent() const
0535    {
0536       return fParent;
0537    }
0538    virtual Double_t GetParameter(Int_t ipar) const
0539    {
0540       return (fFormula) ? fFormula->GetParameter(ipar) : fParams->GetParameter(ipar);
0541    }
0542    virtual Double_t GetParameter(const TString &name)  const
0543    {
0544       return (fFormula) ? fFormula->GetParameter(name) : fParams->GetParameter(name);
0545    }
0546    virtual Double_t *GetParameters() const
0547    {
0548       return (fFormula) ? fFormula->GetParameters() : const_cast<Double_t *>(fParams->GetParameters());
0549    }
0550    virtual void     GetParameters(Double_t *params)
0551    {
0552       if (fFormula) fFormula->GetParameters(params);
0553       else std::copy(fParams->ParamsVec().begin(), fParams->ParamsVec().end(), params);
0554    }
0555    virtual const char *GetParName(Int_t ipar) const
0556    {
0557       return (fFormula) ? fFormula->GetParName(ipar) : fParams->GetParName(ipar);
0558    }
0559    virtual Int_t    GetParNumber(const char *name) const
0560    {
0561       return (fFormula) ? fFormula->GetParNumber(name) : fParams->GetParNumber(name);
0562    }
0563    virtual Double_t GetParError(Int_t ipar) const;
0564    virtual const Double_t *GetParErrors() const
0565    {
0566       return fParErrors.data();
0567    }
0568    virtual void     GetParLimits(Int_t ipar, Double_t &parmin, Double_t &parmax) const;
0569    virtual Double_t GetProb() const;
0570    virtual Int_t    GetQuantiles(Int_t nprobSum, Double_t *q, const Double_t *probSum);
0571    virtual Double_t GetRandom(TRandom * rng = nullptr, Option_t * opt = nullptr);
0572    virtual Double_t GetRandom(Double_t xmin, Double_t xmax, TRandom * rng = nullptr, Option_t * opt = nullptr);
0573    virtual void     GetRange(Double_t &xmin, Double_t &xmax) const;
0574    virtual void     GetRange(Double_t &xmin, Double_t &ymin, Double_t &xmax, Double_t &ymax) const;
0575    virtual void     GetRange(Double_t &xmin, Double_t &ymin, Double_t &zmin, Double_t &xmax, Double_t &ymax, Double_t &zmax) const;
0576    virtual Double_t GetSave(const Double_t *x);
0577    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;
0578    virtual Double_t GetXmin() const
0579    {
0580       return fXmin;
0581    }
0582    virtual Double_t GetXmax() const
0583    {
0584       return fXmax;
0585    }
0586    TAxis           *GetXaxis() const ;
0587    TAxis           *GetYaxis() const ;
0588    TAxis           *GetZaxis() const ;
0589    virtual Double_t GetVariable(const TString &name)
0590    {
0591       return (fFormula) ? fFormula->GetVariable(name) : 0;
0592    }
0593    virtual Double_t GradientPar(Int_t ipar, const Double_t *x, Double_t eps = 0.01);
0594    template <class T>
0595    T GradientPar(Int_t ipar, const T *x, Double_t eps = 0.01);
0596    template <class T>
0597    T GradientParTempl(Int_t ipar, const T *x, Double_t eps = 0.01);
0598
0599    virtual void     GradientPar(const Double_t *x, Double_t *grad, Double_t eps = 0.01);
0600    template <class T>
0601    void GradientPar(const T *x, T *grad, Double_t eps = 0.01);
0602    template <class T>
0603    void GradientParTempl(const T *x, T *grad, Double_t eps = 0.01);
0604
0605    virtual void     InitArgs(const Double_t *x, const Double_t *params);
0606    static  void     InitStandardFunctions();
0607    virtual Double_t Integral(Double_t a, Double_t b, Double_t epsrel = 1.e-12);
0608    virtual Double_t IntegralOneDim(Double_t a, Double_t b, Double_t epsrel, Double_t epsabs, Double_t &err);
0609    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);
0610    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);
0611    // virtual Double_t IntegralFast(const TGraph *g, Double_t a, Double_t b, Double_t *params = nullptr);
0612    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);
0613    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);
0614    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)
0615    {
0616       return  IntegralMultiple(n, a, b, maxpts, epsrel, epsrel, relerr, nfnevl, ifail);
0617    }
0618    virtual Double_t IntegralMultiple(Int_t n, const Double_t *a, const Double_t *b, Double_t epsrel, Double_t &relerr);
0619    virtual Bool_t   IsEvalNormalized() const
0620    {
0621       return fNormalized;
0622    }
0623    /// return kTRUE if the point is inside the function range
0624    virtual Bool_t   IsInside(const Double_t *x) const
0625    {
0626       return !((x[0] < fXmin) || (x[0] > fXmax));
0627    }
0628    virtual Bool_t   IsLinear() const
0629    {
0630       return (fFormula) ? fFormula->IsLinear() : false;
0631    }
0632    virtual Bool_t   IsValid() const;
0633    void     Print(Option_t *option = "") const override;
0634    void     Paint(Option_t *option = "") override;
0635    virtual void     ReleaseParameter(Int_t ipar);
0636    virtual void     Save(Double_t xmin, Double_t xmax, Double_t ymin, Double_t ymax, Double_t zmin, Double_t zmax);
0637    void     SavePrimitive(std::ostream &out, Option_t *option = "") override;
0638    virtual void     SetChisquare(Double_t chi2)
0639    {
0640       fChisquare = chi2;
0641    }
0642    virtual void     SetFitResult(const ROOT::Fit::FitResult &result, const Int_t *indpar = nullptr);
0643    template <class PtrObj, typename MemFn>
0644    void SetFunction(PtrObj &p, MemFn memFn);
0645    template <typename Func>
0646    void SetFunction(Func f);
0647    virtual void     SetMaximum(Double_t maximum = -1111); // *MENU*
0648    virtual void     SetMinimum(Double_t minimum = -1111); // *MENU*
0649    virtual void     SetNDF(Int_t ndf);
0650    virtual void     SetNumberFitPoints(Int_t npfits)
0651    {
0652       fNpfits = npfits;
0653    }
0654    virtual void     SetNormalized(Bool_t flag)
0655    {
0656       fNormalized = flag;
0657       Update();
0658    }
0659    virtual void     SetNpx(Int_t npx = 100); // *MENU*
0660    virtual void     SetParameter(Int_t param, Double_t value)
0661    {
0662       (fFormula) ? fFormula->SetParameter(param, value) : fParams->SetParameter(param, value);
0663       Update();
0664    }
0665    virtual void     SetParameter(const TString &name, Double_t value)
0666    {
0667       (fFormula) ? fFormula->SetParameter(name, value) : fParams->SetParameter(name, value);
0668       Update();
0669    }
0670    virtual void     SetParameters(const Double_t *params)
0671    {
0672       (fFormula) ? fFormula->SetParameters(params) : fParams->SetParameters(params);
0673       Update();
0674    }
0675    /// Set parameter values.
0676    /// NaN values will be skipped, meaning that the corresponding parameters will not be changed.
0677    virtual void SetParameters(double p0, double p1 = TMath::QuietNaN(), double p2 = TMath::QuietNaN(),
0678                               double p3 = TMath::QuietNaN(), double p4 = TMath::QuietNaN(), double p5 = TMath::QuietNaN(),
0679                               double p6 = TMath::QuietNaN(), double p7 = TMath::QuietNaN(), double p8 = TMath::QuietNaN(),
0680                               double p9 = TMath::QuietNaN(), double p10 = TMath::QuietNaN())
0681    {
0682       // Note: this is not made a variadic template method because it would
0683       // presumably break the context menu in the TBrowser. Also, probably this
0684       // method should not be virtual, because if the user wants to change
0685       // parameter setting behavior, the SetParameter() method can be
0686       // overridden.
0687       if (fFormula) fFormula->SetParameters(p0, p1, p2, p3, p4, p5, p6, p7, p8, p9, p10);
0688       else          fParams->SetParameters(p0, p1, p2, p3, p4, p5, p6, p7, p8, p9, p10);
0689       Update();
0690    } // *MENU*
0691    virtual void     SetParName(Int_t ipar, const char *name);
0692    virtual void     SetParNames(const char *name0 = "", const char *name1 = "", const char *name2 = "",
0693                                 const char *name3 = "", const char *name4 = "", const char *name5 = "",
0694                                 const char *name6 = "", const char *name7 = "", const char *name8 = "",
0695                                 const char *name9 = "", const char *name10 = ""); // *MENU*
0696    virtual void     SetParError(Int_t ipar, Double_t error);
0697    virtual void     SetParErrors(const Double_t *errors);
0698    virtual void     SetParLimits(Int_t ipar, Double_t parmin, Double_t parmax);
0699    virtual void     SetParent(TObject *p = nullptr)
0700    {
0701       fParent = p;
0702    }
0703    virtual void     SetRange(Double_t xmin, Double_t xmax); // *MENU*
0704    virtual void     SetRange(Double_t xmin, Double_t ymin,  Double_t xmax, Double_t ymax);
0705    virtual void     SetRange(Double_t xmin, Double_t ymin, Double_t zmin,  Double_t xmax, Double_t ymax, Double_t zmax);
0706    virtual void     SetSavedPoint(Int_t point, Double_t value);
0707    void     SetTitle(const char *title = "") override; // *MENU*
0708    virtual void     SetVectorized(Bool_t vectorized)
0709    {
0710       if (fType == EFType::kFormula && fFormula)
0711          fFormula->SetVectorized(vectorized);
0712       else
0713          Warning("SetVectorized", "Can only set vectorized flag on formula-based TF1");
0714    }
0715    virtual void     Update();
0716
0717    static  TF1     *GetCurrent();
0718    static  void     AbsValue(Bool_t reject = kTRUE);
0719    static  void     RejectPoint(Bool_t reject = kTRUE);
0720    static  Bool_t   RejectedPoint();
0721    static  void     SetCurrent(TF1 *f1);
0722
0723    //Moments
0724    virtual Double_t Moment(Double_t n, Double_t a, Double_t b, const Double_t *params = nullptr, Double_t epsilon = 0.000001);
0725    virtual Double_t CentralMoment(Double_t n, Double_t a, Double_t b, const Double_t *params = nullptr, Double_t epsilon = 0.000001);
0726    virtual Double_t Mean(Double_t a, Double_t b, const Double_t *params = nullptr, Double_t epsilon = 0.000001)
0727    {
0728       return Moment(1, a, b, params, epsilon);
0729    }
0730    virtual Double_t Variance(Double_t a, Double_t b, const Double_t *params = nullptr, Double_t epsilon = 0.000001)
0731    {
0732       return CentralMoment(2, a, b, params, epsilon);
0733    }
0734
0735    //some useful static utility functions to compute sampling points for Integral
0736    //static  void     CalcGaussLegendreSamplingPoints(TGraph *g, Double_t eps=3.0e-11);
0737    //static  TGraph  *CalcGaussLegendreSamplingPoints(Int_t num=21, Double_t eps=3.0e-11);
0738    static  void     CalcGaussLegendreSamplingPoints(Int_t num, Double_t *x, Double_t *w, Double_t eps = 3.0e-11);
0739
0740 private:
0741    template <class T>
0742    T EvalParTempl(const T *data, const Double_t *params = nullptr);
0743
0744 #ifdef R__HAS_VECCORE
0745    inline double EvalParVec(const Double_t *data, const Double_t *params);
0746 #endif
0747
0748    ClassDefOverride(TF1, 12) // The Parametric 1-D function
0749 };
0750
0751 namespace ROOT {
0752    namespace Internal {
0753
0754       template<class Func>
0755       void TF1Builder<Func>::Build(TF1 *f, Func func)
0756       {
0757          using Fnc_t = typename ROOT::Internal::GetFunctorType<decltype(ROOT::Internal::GetTheRightOp(&Func::operator()))>::type;
0758          f->fType = std::is_same<Fnc_t, double>::value? TF1::EFType::kTemplScalar : TF1::EFType::kTemplVec;
0759          f->fFunctor.reset(new TF1::TF1FunctorPointerImpl<Fnc_t>(ROOT::Math::ParamFunctorTempl<Fnc_t>(func)));
0760          f->fParams = std::make_unique<TF1Parameters>(f->fNpar);
0761       }
0762
0763       template<class Func>
0764       void TF1Builder<Func *>::Build(TF1 *f, Func *func)
0765       {
0766          using Fnc_t = typename ROOT::Internal::GetFunctorType<decltype(ROOT::Internal::GetTheRightOp(&Func::operator()))>::type;
0767          f->fType = std::is_same<Fnc_t, double>::value? TF1::EFType::kTemplScalar : TF1::EFType::kTemplVec;
0768          f->fFunctor.reset(new TF1::TF1FunctorPointerImpl<Fnc_t>(ROOT::Math::ParamFunctorTempl<Fnc_t>(func)));
0769          f->fParams = std::make_unique<TF1Parameters>(f->fNpar);
0770       }
0771
0772       /// TF1 building from a string
0773       /// used to build a TFormula based on a lambda function
0774       template<>
0775       struct TF1Builder<const char *> {
0776          static void Build(TF1 *f, const char *formula)
0777          {
0778             f->fType = TF1::EFType::kFormula;
0779             f->fFormula = std::make_unique<TFormula>("tf1lambda", formula, f->fNdim, f->fNpar, false);
0780             TString formulaExpression(formula);
0781             Ssiz_t first = formulaExpression.Index("return") + 7;
0782             Ssiz_t last  = formulaExpression.Last(';');
0783             TString title = formulaExpression(first, last - first);
0784             f->SetTitle(title);
0785          }
0786       };
0787    }
0788 }
0789
0790 inline Double_t TF1::operator()(Double_t x, Double_t y, Double_t z, Double_t t) const
0791 {
0792    return Eval(x, y, z, t);
0793 }
0794
0795 template <class T>
0796 inline T TF1::operator()(const T *x, const Double_t *params)
0797 {
0798    return EvalPar(x, params);
0799 }
0800
0801 ////////////////////////////////////////////////////////////////////////////////
0802 ///   EvalPar for vectorized
0803 template <class T>
0804 T TF1::EvalPar(const T *x, const Double_t *params)
0805 {
0806   if (fType == EFType::kTemplVec || fType == EFType::kTemplScalar) {
0807      return EvalParTempl(x, params);
0808   } else if (fType == EFType::kFormula) {
0809      return fFormula->EvalPar(x, params);
0810   } else
0811      return TF1::EvalPar((double *)x, params);
0812 }
0813
0814 ////////////////////////////////////////////////////////////////////////////////
0815 ///   Eval for vectorized functions
0816 // template <class T>
0817 // T TF1::Eval(T x, T y, T z, T t) const
0818 // {
0819 //    if (fType == EFType::kFormula)
0820 //       return fFormula->Eval(x, y, z, t);
0821
0822 //    T xx[] = {x, y, z, t};
0823 //    Double_t *pp = (Double_t *)fParams->GetParameters();
0824 //    return ((TF1 *)this)->EvalPar(xx, pp);
0825 // }
0826
0827 // Internal to TF1. Evaluates Templated interfaces
0828 template <class T>
0829 inline T TF1::EvalParTempl(const T *data, const Double_t *params)
0830 {
0831    assert(fType == EFType::kTemplScalar || fType == EFType::kTemplVec);
0832    if (!params) params = (Double_t *)fParams->GetParameters();
0833    if (fFunctor)
0834       return ((TF1FunctorPointerImpl<T> *)fFunctor.get())->fImpl(data, params);
0835
0836    // this should throw an error
0837    // we nned to implement a vectorized GetSave(x)
0838    return TMath::SignalingNaN();
0839 }
0840
0841 #ifdef R__HAS_VECCORE
0842 // Internal to TF1. Evaluates Vectorized TF1 on data of type Double_v
0843 inline double TF1::EvalParVec(const Double_t *data, const Double_t *params)
0844 {
0845    assert(fType == EFType::kTemplVec);
0846    std::vector<ROOT::Double_v> d(fNdim);
0847    ROOT::Double_v res;
0848
0849    for(auto i=0; i<fNdim; i++) {
0850       d[i] = ROOT::Double_v(data[i]);
0851    }
0852
0853    if (fFunctor) {
0854       res = ((TF1FunctorPointerImpl<ROOT::Double_v> *) fFunctor.get())->fImpl(d.data(), params);
0855    } else {
0856       //    res = GetSave(x);
0857       return TMath::SignalingNaN();
0858    }
0859    return vecCore::Get<ROOT::Double_v>(res, 0);
0860 }
0861 #endif
0862
0863 inline void TF1::SetRange(Double_t xmin, Double_t,  Double_t xmax, Double_t)
0864 {
0865    TF1::SetRange(xmin, xmax);
0866 }
0867 inline void TF1::SetRange(Double_t xmin, Double_t, Double_t,  Double_t xmax, Double_t, Double_t)
0868 {
0869    TF1::SetRange(xmin, xmax);
0870 }
0871
0872 template <typename Func>
0873 void TF1::SetFunction(Func f)
0874 {
0875    // set function from a generic C++ callable object
0876    fType = EFType::kPtrScalarFreeFcn;
0877    fFunctor = std::make_unique<TF1::TF1FunctorPointerImpl<double>>(ROOT::Math::ParamFunctor(f));
0878 }
0879 template <class PtrObj, typename MemFn>
0880 void TF1::SetFunction(PtrObj &p, MemFn memFn)
0881 {
0882    // set from a pointer to a member function
0883    fType = EFType::kPtrScalarFreeFcn;
0884    fFunctor = std::make_unique<TF1::TF1FunctorPointerImpl<double>>(ROOT::Math::ParamFunctor(p, memFn));
0885 }
0886
0887 template <class T>
0888 inline T TF1::GradientPar(Int_t ipar, const T *x, Double_t eps)
0889 {
0890    if (fType == EFType::kTemplVec || fType == EFType::kTemplScalar) {
0891       return GradientParTempl<T>(ipar, x, eps);
0892    } else
0893       return GradientParTempl<Double_t>(ipar, (const Double_t *)x, eps);
0894 }
0895
0896 template <class T>
0897 inline T TF1::GradientParTempl(Int_t ipar, const T *x, Double_t eps)
0898 {
0899    if (GetNpar() == 0)
0900       return 0;
0901
0902    if (eps < 1e-10 || eps > 1) {
0903       Warning("Derivative", "parameter esp=%g out of allowed range[1e-10,1], reset to 0.01", eps);
0904       eps = 0.01;
0905    }
0906    Double_t h;
0907    TF1 *func = (TF1 *)this;
0908    Double_t *parameters = GetParameters();
0909
0910    // Copy parameters for thread safety
0911    std::vector<Double_t> parametersCopy(parameters, parameters + GetNpar());
0912    parameters = parametersCopy.data();
0913
0914    Double_t al, bl, h2;
0915    T f1, f2, g1, g2, d0, d2;
0916
0917    ((TF1 *)this)->GetParLimits(ipar, al, bl);
0918    if (al * bl != 0 && al >= bl) {
0919       // this parameter is fixed
0920       return 0;
0921    }
0922
0923    // check if error has been computer (is not zero)
0924    if (func->GetParError(ipar) != 0)
0925       h = eps * func->GetParError(ipar);
0926    else
0927       h = eps;
0928
0929    // save original parameters
0930    Double_t par0 = parameters[ipar];
0931
0932    parameters[ipar] = par0 + h;
0933    f1 = func->EvalPar(x, parameters);
0934    parameters[ipar] = par0 - h;
0935    f2 = func->EvalPar(x, parameters);
0936    parameters[ipar] = par0 + h / 2;
0937    g1 = func->EvalPar(x, parameters);
0938    parameters[ipar] = par0 - h / 2;
0939    g2 = func->EvalPar(x, parameters);
0940
0941    // compute the central differences
0942    h2 = 1 / (2. * h);
0943    d0 = f1 - f2;
0944    d2 = 2 * (g1 - g2);
0945
0946    T grad = h2 * (4 * d2 - d0) / 3.;
0947
0948    // restore original value
0949    parameters[ipar] = par0;
0950
0951    return grad;
0952 }
0953
0954 template <class T>
0955 inline void TF1::GradientPar(const T *x, T *grad, Double_t eps)
0956 {
0957    if (fType == EFType::kTemplVec || fType == EFType::kTemplScalar) {
0958       GradientParTempl<T>(x, grad, eps);
0959    } else
0960       GradientParTempl<Double_t>((const Double_t *)x, (Double_t *)grad, eps);
0961 }
0962
0963 template <class T>
0964 inline void TF1::GradientParTempl(const T *x, T *grad, Double_t eps)
0965 {
0966    if (eps < 1e-10 || eps > 1) {
0967       Warning("Derivative", "parameter esp=%g out of allowed range[1e-10,1], reset to 0.01", eps);
0968       eps = 0.01;
0969    }
0970
0971    for (Int_t ipar = 0; ipar < GetNpar(); ipar++) {
0972       grad[ipar] = GradientParTempl<T>(ipar, x, eps);
0973    }
0974 }
0975
0976 #endif