Back to home page

EIC code displayed by LXR
 
 

    

Warning, file /include/root/TF1.h was not indexed or was modified since last indexation (in which case cross-reference links may be missing, inaccurate or erroneous).

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

This page was automatically generated by the 2.3.7 LXR engine.

The LXR team
Valid CSS 2.1! Valid HTML 4.01!