    // @(#)root/hist:$Name:  $:$Id: TFormulaPrimitive.cxx,v 1.2 2005/05/18 12:31:09 brun Exp $
    // Author: Marian Ivanov, 2005

    /*************************************************************************
    * Copyright (C) 1995-2000, Rene Brun and Fons Rademakers.               *
    * All rights reserved.                                                  *
    *                                                                       *
    * For the licensing terms see $ROOTSYS/LICENSE.                         *
    * For the list of contributors see $ROOTSYS/README/CREDITS.             *
    *************************************************************************/

    #include <math.h>

    #include "TFormulaPrimitive.h"
    #include "TMath.h"
    #ifdef WIN32
    #pragma optimize("",off)
    #endif

    void TMath_GenerInterface();

    ClassImp(TFormulaPrimitive)

    //______________________________________________________________________________
    //*-*-*-*-*-*-*-*-*-*-*The  F O R M U L A   P R I M I T I V E   class*-*-*-*-*-*-*-*-*-*-*-*-*-*-*
    //*-*                  =========================
    //*-*
    //*-*
    //*-*    Helper class for TFormula to speed up TFormula evaluation
    //*-*    TFormula can use all functions registered in the list of TFormulaPrimitives
    //*-*    User can add new function to the list of primitives
    //*-*    if FormulaPrimitive with given name is already defined new primitive is ignored
    //*-*    Example:
    //*-*      TFormulaPrimitive::AddFormula(new TFormulaPrimitive("Pow2","Pow2",TFastFun::Pow2));
    //*-*      TF1 f1("f1","Pow2(x)");
    //*-*
    //*-*
    //*-*

    //*-*    TFormulaPrimitive is used to get direct acces to the function pointers
    //*-*    GenFunc     -  pointers  to the static function
    //*-*    TFunc       -  pointers  to the data member functions
    //*-*
    //*-*    The following sufixes are currently used, to describe function arguments:
    //*-*    ------------------------------------------------------------------------
    //*-*    G     - generic layout - pointer to double (arguments), pointer to double (parameters)
    //*-*    10    - double
    //*-*    110   - double, double
    //*-*    1110  - double, double, double
    //*-*

    //*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*

    //_______________________________________________________________________
    // TFormula primitive
    //
    TObjArray * TFormulaPrimitive::fgListOfFunction = 0;

    //______________________________________________________________________________
     TFormulaPrimitive::TFormulaPrimitive() : TNamed()
    {
       //
       //
       fType=0;
       fNArguments=0;
       fNParameters=0;
       fIsStatic=kTRUE;
    }

    //______________________________________________________________________________
     TFormulaPrimitive::TFormulaPrimitive(const char *name,const char *formula,
                                         GenFunc0 fpointer) : TNamed(name,formula)
    {
       fType= 0;
       fFunc0 = fpointer;
       fNArguments=0;
       fNParameters=0;
       fIsStatic=kTRUE;
    }

    //______________________________________________________________________________
     TFormulaPrimitive::TFormulaPrimitive(const char *name,const char *formula,
                                         GenFunc10 fpointer) : TNamed(name,formula)
    {
       //
       //
       //
       fType= 10;
       fFunc10 = fpointer;
       fNArguments=1;
       fNParameters=0;
       fIsStatic=kTRUE;
    }

    //______________________________________________________________________________
     TFormulaPrimitive::TFormulaPrimitive(const char *name,const char *formula,
                                         GenFunc110 fpointer) : TNamed(name,formula)
    {
       //
       //
       //
       fType= 110;
       fFunc110 = fpointer;
       fNArguments=2;
       fNParameters=0;
       fIsStatic=kTRUE;
    }

    //______________________________________________________________________________
     TFormulaPrimitive::TFormulaPrimitive(const char *name,const char *formula,
                                         GenFunc1110 fpointer) : TNamed(name,formula)
    {
       //
       //
       //
       fType= 1110;
       fFunc1110 = fpointer;
       fNArguments=3;
       fNParameters=0;
       fIsStatic=kTRUE;
    }

    //______________________________________________________________________________
     TFormulaPrimitive::TFormulaPrimitive(const char *name,const char *formula,
                                         GenFuncG fpointer,Int_t npar) : TNamed(name,formula)
    {
       //
       //
       //
       fType= -1;
       fFuncG = fpointer;
       fNArguments=2;
       fNParameters=npar;
       fIsStatic=kTRUE;
    }

    //______________________________________________________________________________
     TFormulaPrimitive::TFormulaPrimitive(const char *name,const char *formula,
                                         TFuncG fpointer) : TNamed(name,formula)
    {
       //
       //
       //
       fType= 0;
       fTFuncG=fpointer;
       fNArguments=0;
       fNParameters=0;
       fIsStatic =kFALSE;
    }

    //______________________________________________________________________________
     TFormulaPrimitive::TFormulaPrimitive(const char *name,const char *formula,
                                         TFunc0 fpointer) : TNamed(name,formula)
    {
       //
       //
       //
       fType  = 0;
       fTFunc0= fpointer;
       fNArguments=0;
       fNParameters=0;
       fIsStatic =kFALSE;
    }

    //______________________________________________________________________________
     TFormulaPrimitive::TFormulaPrimitive(const char *name,const char *formula,
                                         TFunc10 fpointer) : TNamed(name,formula)
    {
       //
       //
       //
       fType=-10;
       fTFunc10=fpointer;
       fNArguments=1;
       fNParameters=0;
       fIsStatic =kFALSE;
    }

    //______________________________________________________________________________
     TFormulaPrimitive::TFormulaPrimitive(const char *name,const char *formula,
                                         TFunc110 fpointer) : TNamed(name,formula)
    {
       //
       //
       //
       fType=-110;
       fTFunc110=fpointer;
       fNArguments=2;
       fNParameters=0;
       fIsStatic =kFALSE;
    }

    //______________________________________________________________________________
     TFormulaPrimitive::TFormulaPrimitive(const char *name,const char *formula,
                                         TFunc1110 fpointer) :TNamed(name,formula)
    {
       //
       //
       //
       fType=-1110;
       fTFunc1110=fpointer;
       fNArguments=3;
       fNParameters=0;
       fIsStatic =kFALSE;
    }

    //______________________________________________________________________________
     TFormulaPrimitive* TFormulaPrimitive::FindFormula(const char* name)
    {
       //
       //find the formula in the list of formulas
       //
       if (!fgListOfFunction) {
          BuildBasicFormulas();
       }
       return (TFormulaPrimitive*)fgListOfFunction->FindObject(name);
    }

    //______________________________________________________________________________
     Int_t TFormulaPrimitive::AddFormula(TFormulaPrimitive * formula)
    {
       //
       // add formula to the list of primitive formulas
       // if primitive formula already defined do nothing
       //
       if (fgListOfFunction == 0) BuildBasicFormulas();
       if (FindFormula(formula->GetName())){
          return 0;
       }
       fgListOfFunction->AddLast(formula);
       return 1;
    }

    //______________________________________________________________________________
     Double_t TFormulaPrimitive::Eval(Double_t* x)
    {
       //
       // Eval primitive function at point x
       //
       if (fIsStatic ==kFALSE) return 0;

       if (fType==0) return  fFunc0();
       if (fType==10) {
          return fFunc10(x[0]);
       }
       if (fType==110) {
          return fFunc110(x[0],x[1]);
       }
       if (fType==1110) {
          return fFunc1110(x[0],x[1],x[2]);
       }
       return 0;
    }

    //______________________________________________________________________________
     Double_t  TFormulaPrimitive::Eval(TObject *o, Double_t *x)
    {
       //
       // Eval member function of object o at point x
       //
       if (fIsStatic ==kTRUE) return 0;
       if (fType== 0)    return (*o.*fTFunc0)();
       if (fType==-10)   return (*o.*fTFunc10)(*x);
       if (fType==-110)  return (*o.*fTFunc110)(x[0],x[1]);
       if (fType==-1110) return (*o.*fTFunc1110)(x[0],x[1],x[2]);
       return 0;
    }

    //______________________________________________________________________________
     Double_t TFormulaPrimitive::Eval(Double_t *x, Double_t *param)
    {
       //
       //eval primitive parametric function
       //
       return fFuncG(x,param);
    }

    #define RTFastFun__POLY(var)                                                  \
    {                                                                     \
       Double_t res= param[var-1]+param[var]*x[0];                        \
       for (Int_t j=var-1 ;j>0;j--) res      = param[j-1]+x[0]*res;       \
       return res;                                                        \
    }

    namespace TFastFun {
       //
       // Namespace with basic primitive functions registered by TFormulaPrimitive
       // all function registerd by TFormulaPrimitive can be used in TFormula
       //
       Double_t Pow2(Double_t x){return x*x;}
       Double_t Pow3(Double_t x){return x*x*x;}
       Double_t Pow4(Double_t x){return x*x*x*x;}
       Double_t Pow5(Double_t x){return x*x*x*x*x;}
       inline   Double_t FPoln(Double_t *x, Double_t *param, Int_t npar);
       Double_t FPol0(Double_t * /*x*/, Double_t *param){ return param[0];}
       Double_t FPol1(Double_t *x, Double_t *param){ return param[0]+param[1]*x[0];}
       Double_t FPol2(Double_t *x, Double_t *param){ return param[0]+x[0]*(param[1]+param[2]*x[0]);}
       Double_t FPol3(Double_t *x, Double_t *param){ return param[0]+x[0]*(param[1]+x[0]*(param[2]+param[3]*x[0]));}
       Double_t FPol4(Double_t *x, Double_t *param){ RTFastFun__POLY(4)}
       Double_t FPol5(Double_t *x, Double_t *param){ RTFastFun__POLY(5)}
       Double_t FPol6(Double_t *x, Double_t *param){ RTFastFun__POLY(6)}
       Double_t FPol7(Double_t *x, Double_t *param){ RTFastFun__POLY(7)}
       Double_t FPol8(Double_t *x, Double_t *param){ RTFastFun__POLY(8)}
       Double_t FPol9(Double_t *x, Double_t *param){ RTFastFun__POLY(9)}
       Double_t FPol10(Double_t *x, Double_t *param){ RTFastFun__POLY(10)}
       //
       //
       Double_t PlusXY(Double_t x,Double_t y){return x+y;}
       Double_t MinusXY(Double_t x,Double_t y){return x-y;}
       Double_t MultXY(Double_t x,Double_t y){return x*y;}
       Double_t DivXY(Double_t x, Double_t y){return TMath::Abs(y)>0 ? x/y:0;}
       Double_t XpYpZ(Double_t x, Double_t y, Double_t z){ return x+y+z;}
       Double_t XxYxZ(Double_t x, Double_t y, Double_t z){ return x*y*z;}
       Double_t XxYpZ(Double_t x, Double_t y, Double_t z){ return x*(y+z);}
       Double_t XpYxZ(Double_t x, Double_t y, Double_t z){ return x+(y*z);}
       Double_t Gaus(Double_t x, Double_t mean, Double_t sigma);
       Double_t Gausn(Double_t x, Double_t mean, Double_t sigma);
       Double_t Landau(Double_t x, Double_t mean, Double_t sigma){return TMath::Landau(x,mean,sigma,kFALSE);}
       Double_t Landaun(Double_t x, Double_t mean, Double_t sigma){return TMath::Landau(x,mean,sigma,kTRUE);}
       Double_t Sqrt(Double_t x) {return x>0?sqrt(x):0;}
       //
       Double_t Sign(Double_t x){return (x<0)? -1:1;}
       Double_t Nint(Double_t x){return TMath::Nint(x);}
       Double_t Abs(Double_t x){return TMath::Abs(x);}
       //logical
       Double_t XandY(Double_t x, Double_t y){ return (x*y>0.1);}
       Double_t XorY(Double_t x, Double_t y) { return (x+y>0.1);}
       Double_t XgY(Double_t x, Double_t y) {return (x>y);}
       Double_t XgeY(Double_t x, Double_t y) {return (x>=y);}
       Double_t XlY(Double_t x, Double_t y) {return (x<y);}
       Double_t XleY(Double_t x, Double_t y) {return (x<=y);}
       Double_t XeY(Double_t x,Double_t y) {return (x==y);}
       Double_t XneY(Double_t x,Double_t y) {return (x!=y);}
       Double_t XNot(Double_t x){ return (x<0.1);}
    };

    //______________________________________________________________________________
    Double_t TFastFun::FPoln(Double_t *x, Double_t *param, Int_t npar)
    {
       Double_t res = 0; Double_t temp=1;
       for (Int_t j=npar ;j>=0;j--) {
          res  += temp*param[j];
          temp *= *x;
       }
       return res;
    }

    //______________________________________________________________________________
    Double_t TFastFun::Gaus(Double_t x, Double_t mean, Double_t sigma)
    {
       if (sigma == 0) return 1.e30;
       Double_t arg = (x-mean)/sigma;
       return TMath::Exp(-0.5*arg*arg);
    }

    Double_t TFastFun::Gausn(Double_t x, Double_t mean, Double_t sigma)
    {
       // normalize gauss
       if (sigma == 0)  return 0;
       Double_t arg = (x-mean)/sigma;
       return TMath::Exp(-0.5*arg*arg)/(2.50662827463100024*sigma);  //sqrt(2*Pi)=2.50662827463100024
    }

     Int_t TFormulaPrimitive::BuildBasicFormulas()
    {
       //
       // built-in functions
       //
       if (fgListOfFunction==0) fgListOfFunction = new TObjArray(1000);
       //
       // logical
       //
       AddFormula(new TFormulaPrimitive("XandY","XandY",TFastFun::XandY));
       AddFormula(new TFormulaPrimitive("XorY","XorY",TFastFun::XorY));
       AddFormula(new TFormulaPrimitive("XNot","XNot",TFastFun::XNot));
       AddFormula(new TFormulaPrimitive("XlY","XlY",TFastFun::XlY));
       AddFormula(new TFormulaPrimitive("XleY","XleY",TFastFun::XleY));
       AddFormula(new TFormulaPrimitive("XgY","XgY",TFastFun::XgY));
       AddFormula(new TFormulaPrimitive("XgeY","XgeY",TFastFun::XgeY));
       AddFormula(new TFormulaPrimitive("XeY","XeY",TFastFun::XeY));
       AddFormula(new TFormulaPrimitive("XneY","XneY",TFastFun::XneY));
       // addition  + multiplication
       AddFormula(new TFormulaPrimitive("PlusXY","PlusXY",TFastFun::PlusXY));
       AddFormula(new TFormulaPrimitive("MinusXY","MinusXY",TFastFun::MinusXY));
       AddFormula(new TFormulaPrimitive("MultXY","MultXY",TFastFun::MultXY));
       AddFormula(new TFormulaPrimitive("DivXY","DivXY",TFastFun::DivXY));
       AddFormula(new TFormulaPrimitive("XpYpZ","XpYpZ",TFastFun::XpYpZ));
       AddFormula(new TFormulaPrimitive("XxYxZ","XxYxZ",TFastFun::XxYxZ));
       AddFormula(new TFormulaPrimitive("XxYpZ","XxYpZ",TFastFun::XxYpZ));
       AddFormula(new TFormulaPrimitive("XpYxZ","XpYxZ",TFastFun::XpYxZ));
       //
       //
       AddFormula(new TFormulaPrimitive("Gaus","Gaus",TFastFun::Gaus));
       AddFormula(new TFormulaPrimitive("Gausn","Gausn",TFastFun::Gausn));
       AddFormula(new TFormulaPrimitive("Landau","Landau",TFastFun::Landau));
       AddFormula(new TFormulaPrimitive("Landaun","Landaun",TFastFun::Landau));
       //
       //
       // polynoms
       //
       //
       AddFormula(new TFormulaPrimitive("Pol0","Pol0",(GenFuncG)TFastFun::FPol0,1));
       AddFormula(new TFormulaPrimitive("Pol1","Pol1",(GenFuncG)TFastFun::FPol1,2));
       AddFormula(new TFormulaPrimitive("Pol2","Pol2",(GenFuncG)TFastFun::FPol2,3));
       AddFormula(new TFormulaPrimitive("Pol3","Pol3",(GenFuncG)TFastFun::FPol3,4));
       AddFormula(new TFormulaPrimitive("Pol4","Pol4",(GenFuncG)TFastFun::FPol4,5));
       AddFormula(new TFormulaPrimitive("Pol5","Pol5",(GenFuncG)TFastFun::FPol5,6));
       AddFormula(new TFormulaPrimitive("Pol6","Pol6",(GenFuncG)TFastFun::FPol6,7));
       AddFormula(new TFormulaPrimitive("Pol7","Pol7",(GenFuncG)TFastFun::FPol7,8));
       AddFormula(new TFormulaPrimitive("Pol8","Pol8",(GenFuncG)TFastFun::FPol8,9));
       AddFormula(new TFormulaPrimitive("Pol9","Pol9",(GenFuncG)TFastFun::FPol9,10));
       AddFormula(new TFormulaPrimitive("Pol10","Pol10",(GenFuncG)TFastFun::FPol10,11));
       //
       // pows
       AddFormula(new TFormulaPrimitive("Pow2","Pow2",TFastFun::Pow2));
       AddFormula(new TFormulaPrimitive("Pow3","Pow3",TFastFun::Pow3));
       AddFormula(new TFormulaPrimitive("Pow4","Pow4",TFastFun::Pow4));
       AddFormula(new TFormulaPrimitive("Pow5","Pow5",TFastFun::Pow5));
       //
       //
       AddFormula(new TFormulaPrimitive("TMath::Cos","TMath::Cos",cos));              // 10
       AddFormula(new TFormulaPrimitive("cos","cos",cos));                            // 10
       AddFormula(new TFormulaPrimitive("TMath::Sin","TMath::Sin",sin));              // 11
       AddFormula(new TFormulaPrimitive("sin","sin",sin));                            // 11
       AddFormula(new TFormulaPrimitive("TMath::Tan","TMath::Tan",tan));              // 12
       AddFormula(new TFormulaPrimitive("tan","tan",tan));                            // 12
       AddFormula(new TFormulaPrimitive("TMath::ACos","TMath::ACos",acos));           // 13
       AddFormula(new TFormulaPrimitive("acos","acos",acos));                         // 13
       AddFormula(new TFormulaPrimitive("TMath::ASin","TMath::ASin",asin));           // 14
       AddFormula(new TFormulaPrimitive("asin","asin",asin));                         // 14
       AddFormula(new TFormulaPrimitive("TMath::ATan","TMath::ATan",atan));           // 15
       AddFormula(new TFormulaPrimitive("atan","atan",atan));                         // 15
       AddFormula(new TFormulaPrimitive("TMath::ATan2","TMath::ATan2",atan2));       // 16
       AddFormula(new TFormulaPrimitive("atan2","atan2",atan2));                     // 16
       //   kpow      = 20, ksq = 21, ksqrt     = 22,
       AddFormula(new TFormulaPrimitive("pow","pow",TMath::Power));                         //20
       AddFormula(new TFormulaPrimitive("sq","sq",TFastFun::Pow2));                 //21
       AddFormula(new TFormulaPrimitive("sqrt","sqrt",TFastFun::Sqrt));             //22
       // kmin      = 24, kmax = 25,
       AddFormula(new TFormulaPrimitive("min","min",(GenFunc110)TMath::Min));       //24
       AddFormula(new TFormulaPrimitive("max","max",(GenFunc110)TMath::Max));       //25
       // klog      = 30, kexp = 31, klog10 = 32,
       AddFormula(new TFormulaPrimitive("log","log",TMath::Log));                           //30
       AddFormula(new TFormulaPrimitive("exp","exp",TMath::Exp));                           //31
       AddFormula(new TFormulaPrimitive("log10","log10",TMath::Log10));                     //32
       //
       //    cosh        70                  acosh        73
       //    sinh        71                  asinh        74
       //    tanh        72                  atanh        75
       //
       AddFormula(new TFormulaPrimitive("TMath::CosH","TMath::Cosh",cosh));                     // 70
       AddFormula(new TFormulaPrimitive("cosh","cosh",cosh));                                   // 70
       AddFormula(new TFormulaPrimitive("TMath::SinH","TMath::SinH",sinh));                     // 71
       AddFormula(new TFormulaPrimitive("sinh","sinh",sinh));                                   // 71
       AddFormula(new TFormulaPrimitive("TMath::TanH","TMath::Tanh",tanh));                     // 72
       AddFormula(new TFormulaPrimitive("tanh","tanh",tanh));                                   // 72
       AddFormula(new TFormulaPrimitive("TMath::ACosH","TMath::ACosh",TMath::ACosH));           // 73
       AddFormula(new TFormulaPrimitive("acosh","acosH",TMath::ACosH));                         // 73
       AddFormula(new TFormulaPrimitive("TMath::ASinH","TMath::ASinh",TMath::ASinH));           // 74
       AddFormula(new TFormulaPrimitive("acosh","acosH",TMath::ASinH));                         // 74
       AddFormula(new TFormulaPrimitive("TMath::ATanH","TMath::ATanh",TMath::ATanH));           // 75
       AddFormula(new TFormulaPrimitive("atanh","atanh",TMath::ATanH));                         // 75
       //
       AddFormula(new TFormulaPrimitive("TMath::Abs","TMath::Abs",TMath::Abs));                         //
       AddFormula(new TFormulaPrimitive("TMath::BreitWigner","TMath::BreitWigner",TMath::BreitWigner));                         //
       AddFormula(new TFormulaPrimitive("TMath::Gaus","TMath::Gaus",(GenFunc1110)TMath::Gaus));                         //
       TFormulaPrimitive::AddFormula(new TFormulaPrimitive("TMath::Landau","TMath::Landau",(TFormulaPrimitive::GenFunc1110)TMath::Landau));                         //

       TMath_GenerInterface();
       return 1;
    }