Logo ROOT   6.10/09
Reference Guide
TFormula.cxx
Go to the documentation of this file.
1 // @(#)root/hist:$Id$
2 // Author: Maciej Zimnoch 30/09/2013
3 
4 /*************************************************************************
5  * Copyright (C) 1995-2013, Rene Brun and Fons Rademakers. *
6  * All rights reserved. *
7  * *
8  * For the licensing terms see $ROOTSYS/LICENSE. *
9  * For the list of contributors see $ROOTSYS/README/CREDITS. *
10  *************************************************************************/
11 
12 #if __cplusplus >= 201103L
13 #define ROOT_CPLUSPLUS11 1
14 #endif
15 
16 #include "TROOT.h"
17 #include "TClass.h"
18 #include "TMethod.h"
19 #include "TMath.h"
20 #include "TF1.h"
21 #include "TMethodCall.h"
22 #include <TBenchmark.h>
23 #include "TError.h"
24 #include "TInterpreter.h"
25 #include "TFormula.h"
26 #include <cassert>
27 #include <iostream>
28 #include <unordered_map>
29 #include <functional>
30 
31 using namespace std;
32 
33 // #define __STDC_LIMIT_MACROS
34 // #define __STDC_CONSTANT_MACROS
35 
36 // #include "cling/Interpreter/Interpreter.h"
37 // #include "cling/Interpreter/Value.h"
38 // #include "cling/Interpreter/StoredValueRef.h"
39 
40 
41 #ifdef WIN32
42 #pragma optimize("",off)
43 #endif
44 #include "v5/TFormula.h"
45 
47 
48 /** \class TFormula TFormula.h "inc/TFormula.h"
49  \ingroup Hist
50 The Formula class
51 
52 This is a new version of the TFormula class based on Cling.
53 This class is not 100% backward compatible with the old TFormula class, which is still available in ROOT as
54 `ROOT::v5::TFormula`. Some of the TFormula member funtions available in version 5, such as
55 `Analyze` and `AnalyzeFunction` are not available in the new TFormula.
56 On the other hand formula expressions which were valid in version 5 are still valid in TFormula version 6
57 
58 This class has been implemented during Google Summer of Code 2013 by Maciej Zimnoch.
59 
60 ### Example of valid expressions:
61 
62  - `sin(x)/x`
63  - `[0]*sin(x) + [1]*exp(-[2]*x)`
64  - `x + y**2`
65  - `x^2 + y^2`
66  - `[0]*pow([1],4)`
67  - `2*pi*sqrt(x/y)`
68  - `gaus(0)*expo(3) + ypol3(5)*x`
69  - `gausn(0)*expo(3) + ypol3(5)*x`
70 
71 In the last example above:
72 
73  - `gaus(0)` is a substitute for `[0]*exp(-0.5*((x-[1])/[2])**2)`
74  and (0) means start numbering parameters at 0
75  - `gausn(0)` is a substitute for `[0]*exp(-0.5*((x-[1])/[2])**2)/(sqrt(2*pi)*[2]))`
76  and (0) means start numbering parameters at 0
77  - `expo(3)` is a substitute for `exp([3]+[4]*x`
78  - `pol3(5)` is a substitute for `par[5]+par[6]*x+par[7]*x**2+par[8]*x**3`
79  (`PolN` stands for Polynomial of degree N)
80 
81 `TMath` functions can be part of the expression, eg:
82 
83  - `TMath::Landau(x)*sin(x)`
84  - `TMath::Erf(x)`
85 
86 Formula may contain constans, eg:
87 
88  - `sqrt2`
89  - `e`
90  - `pi`
91  - `ln10`
92  - `infinity`
93 
94 and more.
95 
96 Comparisons operators are also supported `(&amp;&amp;, ||, ==, &lt;=, &gt;=, !)`
97 
98 Examples:
99 
100  `sin(x*(x&lt;0.5 || x&gt;1))`
101 
102 If the result of a comparison is TRUE, the result is 1, otherwise 0.
103 
104 Already predefined names can be given. For example, if the formula
105 
106  `TFormula old("old",sin(x*(x&lt;0.5 || x&gt;1)))`
107 
108 one can assign a name to the formula. By default the name of the object = title = formula itself.
109 
110  `TFormula new("new","x*old")`
111 
112 is equivalent to:
113 
114  `TFormula new("new","x*sin(x*(x&lt;0.5 || x&gt;1))")`
115 
116 The class supports unlimited numer of variables and parameters.
117  By default the names which can be used for the variables are `x,y,z,t` or
118  `x[0],x[1],x[2],x[3],....x[N]` for N-dimensionals formula.
119 
120 This class is not anymore the base class for the function classes `TF1`, but it has now
121 adata member of TF1 which can be access via `TF1::GetFormula`.
122 
123 \class TFormulaFunction
124  Helper class for TFormula
125 
126 \class TFormulaVariable
127  Another helper class for TFormula
128 
129 \class TFormulaParamOrder
130  Functor defining the parameter order
131 */
132 
133 // prefix used for function name passed to Cling
134 static const TString gNamePrefix = "TFormula__";
135 
136 // static map of function pointers and expressions
137 //static std::unordered_map<std::string, TInterpreter::CallFuncIFacePtr_t::Generic_t> gClingFunctions = std::unordered_map<TString, TInterpreter::CallFuncIFacePtr_t::Generic_t>();
138 static std::unordered_map<std::string, void *> gClingFunctions = std::unordered_map<std::string, void * >();
139 
140 ////////////////////////////////////////////////////////////////////////////////
141 Bool_t TFormula::IsOperator(const char c)
142 {
143  // operator ":" must be handled separatly
144  char ops[] = { '+','^', '-','/','*','<','>','|','&','!','=','?'};
145  Int_t opsLen = sizeof(ops)/sizeof(char);
146  for(Int_t i = 0; i < opsLen; ++i)
147  if(ops[i] == c)
148  return true;
149  return false;
150 }
151 
152 ////////////////////////////////////////////////////////////////////////////////
154 {
155  char brackets[] = { ')','(','{','}'};
156  Int_t bracketsLen = sizeof(brackets)/sizeof(char);
157  for(Int_t i = 0; i < bracketsLen; ++i)
158  if(brackets[i] == c)
159  return true;
160  return false;
161 }
162 
163 ////////////////////////////////////////////////////////////////////////////////
165 {
166  return !IsBracket(c) && !IsOperator(c) && c != ',' && c != ' ';
167 }
168 
169 ////////////////////////////////////////////////////////////////////////////////
171 {
172  return name == "x" || name == "z" || name == "y" || name == "t";
173 }
174 
175 ////////////////////////////////////////////////////////////////////////////////
176 Bool_t TFormula::IsScientificNotation(const TString & formula, int i)
177 {
178  // check if the character at position i is part of a scientific notation
179  if ( (formula[i] == 'e' || formula[i] == 'E') && (i > 0 && i < formula.Length()-1) ) {
180  // handle cases: 2e+3 2e-3 2e3 and 2.e+3
181  if ( (isdigit(formula[i-1]) || formula[i-1] == '.') && ( isdigit(formula[i+1]) || formula[i+1] == '+' || formula[i+1] == '-' ) )
182  return true;
183  }
184  return false;
185 }
186 
187 ////////////////////////////////////////////////////////////////////////////////
188 Bool_t TFormula::IsHexadecimal(const TString & formula, int i)
189 {
190  // check if the character at position i is part of a scientific notation
191  if ( (formula[i] == 'x' || formula[i] == 'X') && (i > 0 && i < formula.Length()-1) && formula[i-1] == '0') {
192  if (isdigit(formula[i+1]) )
193  return true;
194  static char hex_values[12] = { 'a','A', 'b','B','c','C','d','D','e','E','f','F'};
195  for (int jjj = 0; jjj < 12; ++jjj) {
196  if (formula[i+1] == hex_values[jjj])
197  return true;
198  }
199  }
200  // else
201  // return false;
202  // // handle cases: 2e+3 2e-3 2e3 and 2.e+3
203  // if ( (isdigit(formula[i-1]) || formula[i-1] == '.') && ( isdigit(formula[i+1]) || formula[i+1] == '+' || formula[i+1] == '-' ) )
204  // return true;
205  // }
206  return false;
207 }
208 ////////////////////////////////////////////////////////////////////////////
209 // check is given position is in a parameter name i.e. within "[ ]"
210 ////
211 Bool_t TFormula::IsAParameterName(const TString & formula, int pos) {
212 
213  Bool_t foundOpenParenthesis = false;
214  if (pos == 0 || pos == formula.Length()-1) return false;
215  for (int i = pos-1; i >=0; i--) {
216  if (formula[i] == ']' ) return false;
217  if (formula[i] == '[' ) {
218  foundOpenParenthesis = true;
219  break;
220  }
221  }
222  if (!foundOpenParenthesis ) return false;
223 
224  // search after the position
225  for (int i = pos+1; i < formula.Length(); i++) {
226  if (formula[i] == ']' ) return true;
227  }
228  return false;
229 }
230 
231 
232 ////////////////////////////////////////////////////////////////////////////////
233 bool TFormulaParamOrder::operator() (const TString& a, const TString& b) const {
234  // implement comparison used to set parameter orders in TFormula
235  // want p2 to be before p10
236 
237  // strip first character in case you have (p0, p1, pN)
238  if ( a[0] == 'p' && a.Length() > 1) {
239  if ( b[0] == 'p' && b.Length() > 1) {
240  // strip first character
241  TString lhs = a(1,a.Length()-1);
242  TString rhs = b(1,b.Length()-1);
243  if (lhs.IsDigit() && rhs.IsDigit() )
244  return (lhs.Atoi() < rhs.Atoi() );
245  }
246  else {
247  return true; // assume a(a numeric name) is always before b (an alphanumeric name)
248  }
249  }
250  else {
251  if ( b[0] == 'p' && b.Length() > 1)
252  // now b is numeric and a is not so return false
253  return false;
254 
255  // case both names are numeric
256  if (a.IsDigit() && b.IsDigit() )
257  return (a.Atoi() < b.Atoi() );
258 
259  }
260 
261  return a < b;
262 }
263 
264 ////////////////////////////////////////////////////////////////////////////////
266 {
267  fName = "";
268  fTitle = "";
269  fClingInput = "";
270  fReadyToExecute = false;
271  fClingInitialized = false;
272  fAllParametersSetted = false;
273  fMethod = 0;
274  fNdim = 0;
275  fNpar = 0;
276  fNumber = 0;
277  fClingName = "";
278  fFormula = "";
279  fLambdaPtr = nullptr;
280 }
281 
282 ////////////////////////////////////////////////////////////////////////////////
283 static bool IsReservedName(const char* name){
284  if (strlen(name)!=1) return false;
285  for (auto const & specialName : {"x","y","z","t"}){
286  if (strcmp(name,specialName)==0) return true;
287  }
288  return false;
289 }
290 
291 ////////////////////////////////////////////////////////////////////////////////
293 {
294 
295  // N.B. a memory leak may happen if user set bit after constructing the object,
296  // Setting of bit should be done only internally
297  if (!TestBit(TFormula::kNotGlobal) && gROOT ) {
299  gROOT->GetListOfFunctions()->Remove(this);
300  }
301 
302  if(fMethod)
303  {
304  fMethod->Delete();
305  }
306  int nLinParts = fLinearParts.size();
307  if (nLinParts > 0) {
308  for (int i = 0; i < nLinParts; ++i) delete fLinearParts[i];
309  }
310 }
311 
312 ////////////////////////////////////////////////////////////////////////////////
313 TFormula::TFormula(const char *name, const char *formula, bool addToGlobList) :
314  TNamed(name,formula),
315  fClingInput(formula),fFormula(formula)
316 {
317  fReadyToExecute = false;
318  fClingInitialized = false;
319  fMethod = 0;
320  fNdim = 0;
321  fNpar = 0;
322  fNumber = 0;
323  fMethod = 0;
324  fLambdaPtr = nullptr;
325 
326  FillDefaults();
327 
328 
329  if (addToGlobList && gROOT) {
330  TFormula *old = 0;
332  old = dynamic_cast<TFormula*> ( gROOT->GetListOfFunctions()->FindObject(name) );
333  if (old)
334  gROOT->GetListOfFunctions()->Remove(old);
335  if (IsReservedName(name))
336  Error("TFormula","The name %s is reserved as a TFormula variable name.\n",name);
337  else
338  gROOT->GetListOfFunctions()->Add(this);
339  }
340  SetBit(kNotGlobal,!addToGlobList);
341 
342  //fName = gNamePrefix + name; // is this needed
343 
344  // do not process null formulas.
345  if (!fFormula.IsNull() ) {
347 
349  }
350 
351 }
352 
353 ////////////////////////////////////////////////////////////////////////////////
354 /// Constructor from a full compileable C++ expression
355 
356 TFormula::TFormula(const char *name, const char *formula, int ndim, int npar, bool addToGlobList) :
357  TNamed(name,formula),
358  fClingInput(formula),fFormula(formula)
359 {
360  fReadyToExecute = false;
361  fClingInitialized = false;
362  fNpar = 0;
363  fMethod = 0;
364  fNumber = 0;
365  fLambdaPtr = nullptr;
366  fFuncPtr = nullptr;
367 
368 
369  fNdim = ndim;
370  for (int i = 0; i < npar; ++i) {
371  DoAddParameter(TString::Format("p%d",i), 0, false);
372  }
373  fAllParametersSetted = true;
374  assert (fNpar == npar);
375 
376  bool ret = InitLambdaExpression(formula);
377 
378  if (ret) {
379 
381 
382  fReadyToExecute = true;
383 
384  if (addToGlobList && gROOT) {
385  TFormula *old = 0;
387  old = dynamic_cast<TFormula*> ( gROOT->GetListOfFunctions()->FindObject(name) );
388  if (old)
389  gROOT->GetListOfFunctions()->Remove(old);
390  if (IsReservedName(name))
391  Error("TFormula","The name %s is reserved as a TFormula variable name.\n",name);
392  else
393  gROOT->GetListOfFunctions()->Add(this);
394  }
395  SetBit(kNotGlobal,!addToGlobList);
396  }
397  else
398  Error("TFormula","Syntax error in building the lambda expression %s", formula );
399 }
400 
401 ////////////////////////////////////////////////////////////////////////////////
402 TFormula::TFormula(const TFormula &formula) :
403  TNamed(formula.GetName(),formula.GetTitle())
404 {
405  fReadyToExecute = false;
406  fClingInitialized = false;
407  fMethod = 0;
408  fNdim = formula.GetNdim();
409  fNpar = formula.GetNpar();
410  fNumber = formula.GetNumber();
411  fFormula = formula.GetExpFormula(); // returns fFormula in case of Lambda's
412  fLambdaPtr = nullptr;
413  fFuncPtr = nullptr;
414 
415  // case of function based on a C++ expression (lambda's) which is ready to be compiled
416  if (formula.fLambdaPtr && formula.TestBit(TFormula::kLambda)) {
417 
419  fParams = formula.fParams;
422 
423  bool ret = InitLambdaExpression(fFormula);
424 
425  if (ret) {
427  fReadyToExecute = true;
428  }
429  else
430  Error("TFormula","Syntax error in building the lambda expression %s", fFormula.Data() );
431 
432  }
433  else {
434 
435  FillDefaults();
436 
439  }
440 
441 
442  if (!TestBit(TFormula::kNotGlobal) && gROOT ) {
444  TFormula *old = (TFormula*)gROOT->GetListOfFunctions()->FindObject(formula.GetName());
445  if (old)
446  gROOT->GetListOfFunctions()->Remove(old);
447 
448  if (IsReservedName(formula.GetName())) {
449  Error("TFormula","The name %s is reserved as a TFormula variable name.\n",formula.GetName());
450  } else
451  gROOT->GetListOfFunctions()->Add(this);
452  }
453 
454 }
455 
456 ////////////////////////////////////////////////////////////////////////////////
457 /// = operator.
458 
460 {
461 
462  if (this != &rhs) {
463  rhs.Copy(*this);
464  }
465  return *this;
466 }
467 
468 ////////////////////////////////////////////////////////////////////////////////
469 Bool_t TFormula::InitLambdaExpression(const char * formula) {
470 
471  std::string lambdaExpression = formula;
472 
473  // check if formula exist already in the map
474  {
476 
477  auto funcit = gClingFunctions.find(lambdaExpression);
478  if (funcit != gClingFunctions.end() ) {
479  fLambdaPtr = funcit->second;
480  fClingInitialized = true;
481  return true;
482  }
483  }
484 
485  // set the cling name using hash of the static formulae map
486  auto hasher = gClingFunctions.hash_function();
487  TString lambdaName = TString::Format("lambda__id%zu", hasher(lambdaExpression) );
488 
489  //lambdaExpression = TString::Format("[&](double * x, double *){ return %s ;}",formula);
490  //TString lambdaName = TString::Format("mylambda_%s",GetName() );
491  TString lineExpr = TString::Format("std::function<double(double*,double*)> %s = %s ;",lambdaName.Data(), lambdaExpression.c_str() );
492  gInterpreter->ProcessLine(lineExpr);
493  fLambdaPtr = (void*) gInterpreter->ProcessLine(TString(lambdaName)+TString(";")); // add ; to avoid printing
494  if (fLambdaPtr != nullptr) {
496  gClingFunctions.insert ( std::make_pair ( lambdaExpression, fLambdaPtr) );
497  fClingInitialized = true;
498  return true;
499  }
500  fClingInitialized = false;
501  return false;
502 }
503 
504 ////////////////////////////////////////////////////////////////////////////////
505 /// Compile the given expression with Cling
506 /// backward compatibility method to be used in combination with the empty constructor
507 /// if no expression is given , the current stored formula (retrieved with GetExpFormula()) or the title is used.
508 /// return 0 if the formula compilation is successfull
509 
510 Int_t TFormula::Compile(const char *expression)
511 {
512  TString formula = expression;
513  if (formula.IsNull() ) {
514  formula = fFormula;
515  if (formula.IsNull() ) formula = GetTitle();
516  }
517 
518  if (formula.IsNull() ) return -1;
519 
520  // do not re-process if it was done before
521  if (IsValid() && formula == fFormula ) return 0;
522 
523  // clear if a formula was already existing
524  if (!fFormula.IsNull() ) Clear();
525 
526  fFormula = formula;
527 
528  if (TestBit(TFormula::kLambda) ) {
529  bool ret = InitLambdaExpression(fFormula);
530  return (ret) ? 0 : 1;
531  }
532 
533  if (fVars.empty() ) FillDefaults();
534  // prepare the formula for Cling
535  //printf("compile: processing formula %s\n",fFormula.Data() );
537  // pass formula in CLing
538  bool ret = PrepareFormula(fFormula);
539 
540  return (ret) ? 0 : 1;
541 }
542 
543 ////////////////////////////////////////////////////////////////////////////////
544 void TFormula::Copy(TObject &obj) const
545 {
546  TNamed::Copy(obj);
547  // need to copy also cling parameters
548  TFormula & fnew = dynamic_cast<TFormula&>(obj);
549 
552 
553  fnew.fFuncs = fFuncs;
554  fnew.fVars = fVars;
555  fnew.fParams = fParams;
556  fnew.fConsts = fConsts;
558  fnew.fFormula = fFormula;
559  fnew.fNdim = fNdim;
560  fnew.fNpar = fNpar;
561  fnew.fNumber = fNumber;
563  // copy Linear parts (it is a vector of TFormula pointers) needs to be copied one by one
564  // looping at all the elements
565  // delete first previous elements
566  int nLinParts = fnew.fLinearParts.size();
567  if (nLinParts > 0) {
568  for (int i = 0; i < nLinParts; ++i) delete fnew.fLinearParts[i];
569  fnew.fLinearParts.clear();
570  }
571  // old size that needs to be copied
572  nLinParts = fLinearParts.size();
573  if (nLinParts > 0) {
574  fnew.fLinearParts.reserve(nLinParts);
575  for (int i = 0; i < nLinParts; ++i) {
576  TFormula * linearNew = new TFormula();
577  TFormula * linearOld = (TFormula*) fLinearParts[i];
578  if (linearOld) {
579  linearOld->Copy(*linearNew);
580  fnew.fLinearParts.push_back(linearNew);
581  }
582  else
583  Warning("Copy","Function %s - expr %s has a dummy linear part %d",GetName(),GetExpFormula().Data(),i);
584  }
585  }
586 
587  fnew.fClingInput = fClingInput;
591  fnew.fClingName = fClingName;
592 
593  // case of function based on a C++ expression (lambda's) which is ready to be compiled
595 
596  bool ret = fnew.InitLambdaExpression(fnew.fFormula);
597  if (ret) {
599  fnew.fReadyToExecute = true;
600  }
601  else {
602  Error("TFormula","Syntax error in building the lambda expression %s", fFormula.Data() );
603  fnew.fReadyToExecute = false;
604  }
605  }
606  else if (fMethod) {
607  if (fnew.fMethod) delete fnew.fMethod;
608  // use copy-constructor of TMethodCall
609  TMethodCall *m = new TMethodCall(*fMethod);
610  fnew.fMethod = m;
611  }
612 
613  fnew.fFuncPtr = fFuncPtr;
614 
615 }
616 
617 ////////////////////////////////////////////////////////////////////////////////
618 /// Clear the formula setting expression to empty and reset the variables and
619 /// parameters containers.
620 
622 {
623  fNdim = 0;
624  fNpar = 0;
625  fNumber = 0;
626  fFormula = "";
627  fClingName = "";
628 
629 
630  if(fMethod) fMethod->Delete();
631  fMethod = nullptr;
632 
633  fClingVariables.clear();
634  fClingParameters.clear();
635  fReadyToExecute = false;
636  fClingInitialized = false;
637  fAllParametersSetted = false;
638  fFuncs.clear();
639  fVars.clear();
640  fParams.clear();
641  fConsts.clear();
642  fFunctionsShortcuts.clear();
643 
644  // delete linear parts
645  int nLinParts = fLinearParts.size();
646  if (nLinParts > 0) {
647  for (int i = 0; i < nLinParts; ++i) delete fLinearParts[i];
648  }
649  fLinearParts.clear();
650 
651 }
652 
653 ////////////////////////////////////////////////////////////////////////////////
654 /// Sets TMethodCall to function inside Cling environment.
655 /// TFormula uses it to execute function.
656 /// After call, TFormula should be ready to evaluate formula.
657 
659 {
660 
661  if(!fMethod)
662  {
663  fMethod = new TMethodCall();
664 
665  Bool_t hasParameters = (fNpar > 0);
666  Bool_t hasVariables = (fNdim > 0);
667  TString prototypeArguments = "";
668  if(hasVariables)
669  {
670  prototypeArguments.Append("Double_t*");
671  }
672  if(hasVariables && hasParameters)
673  {
674  prototypeArguments.Append(",");
675  }
676  if(hasParameters)
677  {
678  prototypeArguments.Append("Double_t*");
679  }
680  // init method call using real function name (cling name) which is defined in ProcessFormula
681  fMethod->InitWithPrototype(fClingName,prototypeArguments);
682  if(!fMethod->IsValid())
683  {
684  Error("Eval","Can't find %s function prototype with arguments %s",fClingName.Data(),prototypeArguments.Data());
685  return false;
686  }
687 
688  // not needed anymore since we use the function pointer
689  // if(hasParameters)
690  // {
691  // Long_t args[2];
692  // args[0] = (Long_t)fClingVariables.data();
693  // args[1] = (Long_t)fClingParameters.data();
694  // fMethod->SetParamPtrs(args,2);
695  // }
696  // else
697  // {
698  // Long_t args[1];
699  // args[0] = (Long_t)fClingVariables.data();
700  // fMethod->SetParamPtrs(args,1);
701  // }
702 
703  CallFunc_t * callfunc = fMethod->GetCallFunc();
705  fFuncPtr = faceptr.fGeneric;
706  }
707  return true;
708 }
709 
710 ////////////////////////////////////////////////////////////////////////////////
711 /// Inputs formula, transfered to C++ code into Cling
712 
714 {
715 
716  if(!fClingInitialized && fReadyToExecute && fClingInput.Length() > 0)
717  {
720  }
721 }
722 
723 ////////////////////////////////////////////////////////////////////////////////
724 /// Fill structures with default variables, constants and function shortcuts
725 
727 {
728 //#ifdef ROOT_CPLUSPLUS11
729 
730  const TString defvars[] = { "x","y","z","t"};
731  const pair<TString,Double_t> defconsts[] = { {"pi",TMath::Pi()}, {"sqrt2",TMath::Sqrt2()},
732  {"infinity",TMath::Infinity()}, {"e",TMath::E()}, {"ln10",TMath::Ln10()},
733  {"loge",TMath::LogE()}, {"c",TMath::C()}, {"g",TMath::G()},
734  {"h",TMath::H()}, {"k",TMath::K()},{"sigma",TMath::Sigma()},
735  {"r",TMath::R()}, {"eg",TMath::EulerGamma()},{"true",1},{"false",0} };
736  // const pair<TString,Double_t> defconsts[] = { {"pi",TMath::Pi()}, {"sqrt2",TMath::Sqrt2()},
737  // {"infinity",TMath::Infinity()}, {"ln10",TMath::Ln10()},
738  // {"loge",TMath::LogE()}, {"true",1},{"false",0} };
739  const pair<TString,TString> funShortcuts[] =
740  { {"sin","TMath::Sin" },
741  {"cos","TMath::Cos" }, {"exp","TMath::Exp"}, {"log","TMath::Log"}, {"log10","TMath::Log10"},
742  {"tan","TMath::Tan"}, {"sinh","TMath::SinH"}, {"cosh","TMath::CosH"},
743  {"tanh","TMath::TanH"}, {"asin","TMath::ASin"}, {"acos","TMath::ACos"},
744  {"atan","TMath::ATan"}, {"atan2","TMath::ATan2"}, {"sqrt","TMath::Sqrt"},
745  {"ceil","TMath::Ceil"}, {"floor","TMath::Floor"}, {"pow","TMath::Power"},
746  {"binomial","TMath::Binomial"},{"abs","TMath::Abs"},
747  {"min","TMath::Min"},{"max","TMath::Max"},{"sign","TMath::Sign" },
748  {"sq","TMath::Sq"}
749  };
750 
751  std::vector<TString> defvars2(10);
752  for (int i = 0; i < 9; ++i)
753  defvars2[i] = TString::Format("x[%d]",i);
754 
755  for(auto var : defvars)
756  {
757  int pos = fVars.size();
758  fVars[var] = TFormulaVariable(var,0,pos);
759  fClingVariables.push_back(0);
760  }
761  // add also the variables definesd like x[0],x[1],x[2],...
762  // support up to x[9] - if needed extend that to higher value
763  // const int maxdim = 10;
764  // for (int i = 0; i < maxdim; ++i) {
765  // TString xvar = TString::Format("x[%d]",i);
766  // fVars[xvar] = TFormulaVariable(xvar,0,i);
767  // fClingVariables.push_back(0);
768  // }
769 
770  for(auto con : defconsts)
771  {
772  fConsts[con.first] = con.second;
773  }
774  for(auto fun : funShortcuts)
775  {
776  fFunctionsShortcuts[fun.first] = fun.second;
777  }
778 
779 /*** - old code tu support C++03
780 #else
781 
782  TString defvarsNames[] = {"x","y","z","t"};
783  Int_t defvarsLength = sizeof(defvarsNames)/sizeof(TString);
784 
785  TString defconstsNames[] = {"pi","sqrt2","infinity","e","ln10","loge","c","g","h","k","sigma","r","eg","true","false"};
786  Double_t defconstsValues[] = {TMath::Pi(),TMath::Sqrt2(),TMath::Infinity(),TMath::E(),TMath::Ln10(),TMath::LogE(),
787  TMath::C(),TMath::G(),TMath::H(),TMath::K(),TMath::Sigma(),TMath::R(),TMath::EulerGamma(), 1, 0};
788  Int_t defconstsLength = sizeof(defconstsNames)/sizeof(TString);
789 
790  TString funShortcutsNames[] = {"sin","cos","exp","log","tan","sinh","cosh","tanh","asin","acos","atan","atan2","sqrt",
791  "ceil","floor","pow","binomial","abs"};
792  TString funShortcutsExtendedNames[] = {"TMath::Sin","TMath::Cos","TMath::Exp","TMath::Log","TMath::Tan","TMath::SinH",
793  "TMath::CosH","TMath::TanH","TMath::ASin","TMath::ACos","TMath::ATan","TMath::ATan2",
794  "TMath::Sqrt","TMath::Ceil","TMath::Floor","TMath::Power","TMath::Binomial","TMath::Abs"};
795  Int_t funShortcutsLength = sizeof(funShortcutsNames)/sizeof(TString);
796 
797  for(Int_t i = 0; i < defvarsLength; ++i)
798  {
799  TString var = defvarsNames[i];
800  Double_t value = 0;
801  unsigned int pos = fVars.size();
802  fVars[var] = TFormulaVariable(var,value,pos);
803  fClingVariables.push_back(value);
804  }
805 
806  for(Int_t i = 0; i < defconstsLength; ++i)
807  {
808  fConsts[defconstsNames[i]] = defconstsValues[i];
809  }
810  for(Int_t i = 0; i < funShortcutsLength; ++i)
811  {
812  pair<TString,TString> fun(funShortcutsNames[i],funShortcutsExtendedNames[i]);
813  fFunctionsShortcuts[fun.first] = fun.second;
814  }
815 
816 #endif
817 ***/
818 
819 }
820 
821 ////////////////////////////////////////////////////////////////////////////////
822 /// Handling polN
823 /// If before 'pol' exist any name, this name will be treated as variable used in polynomial
824 /// eg.
825 /// varpol2(5) will be replaced with: [5] + [6]*var + [7]*var^2
826 /// Empty name is treated like variable x.
827 
828 void TFormula::HandlePolN(TString &formula)
829 {
830  Int_t polPos = formula.Index("pol");
831  while(polPos != kNPOS && !IsAParameterName(formula,polPos) )
832  {
833 
834  Bool_t defaultVariable = false;
835  TString variable;
836  Int_t openingBracketPos = formula.Index('(',polPos);
837  Bool_t defaultCounter = openingBracketPos == kNPOS;
838  Bool_t defaultDegree = true;
839  Int_t degree,counter;
840  TString sdegree;
841  if(!defaultCounter)
842  {
843  // veryfy first of opening parenthesis belongs to pol expression
844  // character between 'pol' and '(' must all be digits
845  sdegree = formula(polPos + 3,openingBracketPos - polPos - 3);
846  if (!sdegree.IsDigit() ) defaultCounter = true;
847  }
848  if (!defaultCounter) {
849  degree = sdegree.Atoi();
850  counter = TString(formula(openingBracketPos+1,formula.Index(')',polPos) - openingBracketPos)).Atoi();
851  }
852  else
853  {
854  Int_t temp = polPos+3;
855  while(temp < formula.Length() && isdigit(formula[temp]))
856  {
857  defaultDegree = false;
858  temp++;
859  }
860  degree = TString(formula(polPos+3,temp - polPos - 3)).Atoi();
861  counter = 0;
862  }
863 
864  TString replacement = TString::Format("[%d]",counter);
865  if(polPos - 1 < 0 || !IsFunctionNameChar(formula[polPos-1]) || formula[polPos-1] == ':' )
866  {
867  variable = "x";
868  defaultVariable = true;
869  }
870  else
871  {
872  Int_t tmp = polPos - 1;
873  while(tmp >= 0 && IsFunctionNameChar(formula[tmp]) && formula[tmp] != ':')
874  {
875  tmp--;
876  }
877  variable = formula(tmp + 1, polPos - (tmp+1));
878  }
879  Int_t param = counter + 1;
880  Int_t tmp = 1;
881  while(tmp <= degree)
882  {
883  if (tmp > 1)
884  replacement.Append(TString::Format("+[%d]*%s^%d",param,variable.Data(),tmp));
885  else
886  replacement.Append(TString::Format("+[%d]*%s",param,variable.Data()));
887  param++;
888  tmp++;
889  }
890  // add paranthesis before and after
891  if (degree > 0) {
892  replacement.Insert(0,'(');
893  replacement.Append(')');
894  }
895  TString pattern;
896  if(defaultCounter && !defaultDegree)
897  {
898  pattern = TString::Format("%spol%d",(defaultVariable ? "" : variable.Data()),degree);
899  }
900  else if(defaultCounter && defaultDegree)
901  {
902  pattern = TString::Format("%spol",(defaultVariable ? "" : variable.Data()));
903  }
904  else
905  {
906  pattern = TString::Format("%spol%d(%d)",(defaultVariable ? "" : variable.Data()),degree,counter);
907  }
908 
909  if (!formula.Contains(pattern)) {
910  Error("HandlePolN","Error handling polynomial function - expression is %s - trying to replace %s with %s ", formula.Data(), pattern.Data(), replacement.Data() );
911  break;
912  }
913  if (formula == pattern) {
914  // case of single polynomial
915  SetBit(kLinear,1);
916  fNumber = 300 + degree;
917  }
918  formula.ReplaceAll(pattern,replacement);
919  polPos = formula.Index("pol");
920  }
921 }
922 
923 ////////////////////////////////////////////////////////////////////////////////
924 /// Handling parametrized functions
925 /// Function can be normalized, and have different variable then x.
926 /// Variables should be placed in brackets after function name.
927 /// No brackets are treated like [x].
928 /// Normalized function has char 'n' after name, eg.
929 /// gausn[var](0) will be replaced with [0]*exp(-0.5*((var-[1])/[2])^2)/(sqrt(2*pi)*[2])
930 ///
931 /// Adding function is easy, just follow these rules:
932 /// - Key for function map is pair of name and dimension of function
933 /// - value of key is a pair function body and normalized function body
934 /// - {Vn} is a place where variable appear, n represents n-th variable from variable list.
935 /// Count starts from 0.
936 /// - [num] stands for parameter number.
937 /// If user pass to function argument 5, num will stand for (5 + num) parameter.
938 ///
939 
941 {
942  map< pair<TString,Int_t> ,pair<TString,TString> > functions;
943  functions.insert(make_pair(make_pair("gaus",1),make_pair("[0]*exp(-0.5*(({V0}-[1])/[2])*(({V0}-[1])/[2]))","[0]*exp(-0.5*(({V0}-[1])/[2])*(({V0}-[1])/[2]))/(sqrt(2*pi)*[2])")));
944  functions.insert(make_pair(make_pair("landau",1),make_pair("[0]*TMath::Landau({V0},[1],[2],false)","[0]*TMath::Landau({V0},[1],[2],true)")));
945  functions.insert(make_pair(make_pair("expo",1),make_pair("exp([0]+[1]*{V0})","")));
946  functions.insert(make_pair(make_pair("crystalball",1),make_pair("[0]*ROOT::Math::crystalball_function({V0},[3],[4],[2],[1])","[0]*ROOT::Math::crystalball_pdf({V0},[3],[4],[2],[1])")));
947  functions.insert(make_pair(make_pair("breitwigner",1),make_pair("[0]*ROOT::Math::breitwigner_pdf({V0},[2],[1])","[0]*ROOT::Math::breitwigner_pdf({V0},[2],[4],[1])")));
948  // chebyshev polynomial
949  functions.insert(make_pair(make_pair("cheb0" ,1),make_pair("ROOT::Math::Chebyshev0({V0},[0])","")));
950  functions.insert(make_pair(make_pair("cheb1" ,1),make_pair("ROOT::Math::Chebyshev1({V0},[0],[1])","")));
951  functions.insert(make_pair(make_pair("cheb2" ,1),make_pair("ROOT::Math::Chebyshev2({V0},[0],[1],[2])","")));
952  functions.insert(make_pair(make_pair("cheb3" ,1),make_pair("ROOT::Math::Chebyshev3({V0},[0],[1],[2],[3])","")));
953  functions.insert(make_pair(make_pair("cheb4" ,1),make_pair("ROOT::Math::Chebyshev4({V0},[0],[1],[2],[3],[4])","")));
954  functions.insert(make_pair(make_pair("cheb5" ,1),make_pair("ROOT::Math::Chebyshev5({V0},[0],[1],[2],[3],[4],[5])","")));
955  functions.insert(make_pair(make_pair("cheb6" ,1),make_pair("ROOT::Math::Chebyshev6({V0},[0],[1],[2],[3],[4],[5],[6])","")));
956  functions.insert(make_pair(make_pair("cheb7" ,1),make_pair("ROOT::Math::Chebyshev7({V0},[0],[1],[2],[3],[4],[5],[6],[7])","")));
957  functions.insert(make_pair(make_pair("cheb8" ,1),make_pair("ROOT::Math::Chebyshev8({V0},[0],[1],[2],[3],[4],[5],[6],[7],[8])","")));
958  functions.insert(make_pair(make_pair("cheb9" ,1),make_pair("ROOT::Math::Chebyshev9({V0},[0],[1],[2],[3],[4],[5],[6],[7],[8],[9])","")));
959  functions.insert(make_pair(make_pair("cheb10",1),make_pair("ROOT::Math::Chebyshev10({V0},[0],[1],[2],[3],[4],[5],[6],[7],[8],[9],[10])","")));
960  // 2-dimensional functions
961  functions.insert(make_pair(make_pair("gaus",2),make_pair("[0]*exp(-0.5*(({V0}-[1])/[2])^2 - 0.5*(({V1}-[3])/[4])^2)","")));
962  functions.insert(make_pair(make_pair("landau",2),make_pair("[0]*TMath::Landau({V0},[1],[2],false)*TMath::Landau({V1},[3],[4],false)","")));
963  functions.insert(make_pair(make_pair("expo",2),make_pair("exp([0]+[1]*{V0})","exp([0]+[1]*{V0}+[2]*{V1})")));
964  // gaussian with correlations
965  functions.insert(make_pair(make_pair("bigaus",2),make_pair("[0]*ROOT::Math::bigaussian_pdf({V0},{V1},[2],[4],[5],[1],[3])","[0]*ROOT::Math::bigaussian_pdf({V0},{V1},[2],[4],[5],[1],[3])")));
966 
967  map<TString,Int_t> functionsNumbers;
968  functionsNumbers["gaus"] = 100;
969  functionsNumbers["bigaus"] = 102;
970  functionsNumbers["landau"] = 400;
971  functionsNumbers["expo"] = 200;
972  functionsNumbers["crystalball"] = 500;
973 
974  // replace old names xygaus -> gaus[x,y]
975  formula.ReplaceAll("xygaus","gaus[x,y]");
976  formula.ReplaceAll("xgaus","gaus[x]");
977  formula.ReplaceAll("ygaus","gaus[y]");
978  formula.ReplaceAll("zgaus","gaus[z]");
979  formula.ReplaceAll("xexpo","expo[x]");
980  formula.ReplaceAll("yexpo","expo[y]");
981  formula.ReplaceAll("zexpo","expo[z]");
982  formula.ReplaceAll("xylandau","landau[x,y]");
983  formula.ReplaceAll("xyexpo","expo[x,y]");
984  // at the moment pre-defined functions have no more than 3 dimensions
985  const char * defaultVariableNames[] = { "x","y","z"};
986 
987  for(map<pair<TString,Int_t>,pair<TString,TString> >::iterator it = functions.begin(); it != functions.end(); ++it)
988  {
989 
990  TString funName = it->first.first;
991  Int_t funDim = it->first.second;
992  Int_t funPos = formula.Index(funName);
993 
994 
995 
996  //std::cout << formula << " ---- " << funName << " " << funPos << std::endl;
997  while(funPos != kNPOS && !IsAParameterName(formula,funPos))
998  {
999 
1000  // should also check that function is not something else (e.g. exponential - parse the expo)
1001  Int_t lastFunPos = funPos + funName.Length();
1002 
1003  // check that first and last character is not alphanumeric
1004  Int_t iposBefore = funPos - 1;
1005  //std::cout << "looping on funpos is " << funPos << " formula is " << formula << " function " << funName << std::endl;
1006  if (iposBefore >= 0) {
1007  assert( iposBefore < formula.Length() );
1008  if (isalpha(formula[iposBefore] ) ) {
1009  //std::cout << "previous character for function " << funName << " is " << formula[iposBefore] << "- skip " << std::endl;
1010  funPos = formula.Index(funName,lastFunPos);
1011  continue;
1012  }
1013  }
1014 
1015  Bool_t isNormalized = false;
1016  if (lastFunPos < formula.Length() ) {
1017  // check if function is normalized by looking at "n" character after function name (e.g. gausn)
1018  isNormalized = (formula[lastFunPos] == 'n');
1019  if (isNormalized) lastFunPos += 1;
1020  if (lastFunPos < formula.Length() ) {
1021  char c = formula[lastFunPos];
1022  // check if also last character is not alphanumeric or is not an operator and not a parenthesis ( or [.
1023  // Parenthesis [] are used to express the variables
1024  if ( isalnum(c ) || ( ! IsOperator(c ) && c != '(' && c != ')' && c != '[' && c != ']' ) ) {
1025  //std::cout << "last character for function " << funName << " is " << c << " skip .." << std::endl;
1026  funPos = formula.Index(funName,lastFunPos);
1027  continue;
1028  }
1029  }
1030  }
1031 
1032  if(isNormalized)
1033  {
1034  SetBit(kNormalized,1);
1035  }
1036  std::vector<TString> variables;
1037  Int_t dim = 0;
1038  TString varList = "";
1039  Bool_t defaultVariables = false;
1040 
1041  // check if function has specified the [...] e.g. gaus[x,y]
1042  Int_t openingBracketPos = funPos + funName.Length() + (isNormalized ? 1 : 0);
1043  Int_t closingBracketPos = kNPOS;
1044  if(openingBracketPos > formula.Length() || formula[openingBracketPos] != '[')
1045  {
1046  dim = funDim;
1047  variables.resize(dim);
1048  for (Int_t idim = 0; idim < dim; ++idim)
1049  variables[idim] = defaultVariableNames[idim];
1050  defaultVariables = true;
1051  }
1052  else
1053  {
1054  // in case of [..] found, assume they specify all the variables. Use it to get function dimension
1055  closingBracketPos = formula.Index(']',openingBracketPos);
1056  varList = formula(openingBracketPos+1,closingBracketPos - openingBracketPos - 1);
1057  dim = varList.CountChar(',') + 1;
1058  variables.resize(dim);
1059  Int_t Nvar = 0;
1060  TString varName = "";
1061  for(Int_t i = 0 ; i < varList.Length(); ++i)
1062  {
1063  if(IsFunctionNameChar(varList[i]))
1064  {
1065  varName.Append(varList[i]);
1066  }
1067  if(varList[i] == ',')
1068  {
1069  variables[Nvar] = varName;
1070  varName = "";
1071  Nvar++;
1072  }
1073  }
1074  if(varName != "") // we will miss last variable
1075  {
1076  variables[Nvar] = varName;
1077  }
1078  }
1079  // chech if dimension obtained from [...] is compatible with what is defined in existing pre-defined functions
1080  //std::cout << " Found dim = " << dim << " and function dimension is " << funDim << std::endl;
1081  if(dim != funDim)
1082  {
1083  pair<TString,Int_t> key = make_pair(funName,dim);
1084  if(functions.find(key) == functions.end())
1085  {
1086  Error("PreProcessFormula","Dimension of function %s is detected to be of dimension %d and is not compatible with existing pre-defined function which has dim %d",
1087  funName.Data(),dim,funDim);
1088  return;
1089  }
1090  // skip the particular function found - we might find later on the corresponding pre-defined function
1091  funPos = formula.Index(funName,lastFunPos);
1092  continue;
1093  }
1094  // look now for the (..) brackets to get the parameter counter (e.g. gaus(0) + gaus(3) )
1095  // need to start for a position
1096  Int_t openingParenthesisPos = (closingBracketPos == kNPOS) ? openingBracketPos : closingBracketPos + 1;
1097  bool defaultCounter = (openingParenthesisPos > formula.Length() || formula[openingParenthesisPos] != '(');
1098 
1099  //Int_t openingParenthesisPos = formula.Index('(',funPos);
1100  //Bool_t defaultCounter = (openingParenthesisPos == kNPOS);
1101  Int_t counter;
1102  if(defaultCounter)
1103  {
1104  counter = 0;
1105  }
1106  else
1107  {
1108  counter = TString(formula(openingParenthesisPos+1,formula.Index(')',funPos) - openingParenthesisPos -1)).Atoi();
1109  }
1110  //std::cout << "openingParenthesisPos " << openingParenthesisPos << " counter is " << counter << std::endl;
1111 
1112  TString body = (isNormalized ? it->second.second : it->second.first);
1113  if(isNormalized && body == "")
1114  {
1115  Error("PreprocessFormula","%d dimension function %s has no normalized form.",it->first.second,funName.Data());
1116  break;
1117  }
1118  for(int i = 0 ; i < body.Length() ; ++i)
1119  {
1120  if(body[i] == '{')
1121  {
1122  // replace {Vn} with variable names
1123  i += 2; // skip '{' and 'V'
1124  Int_t num = TString(body(i,body.Index('}',i) - i)).Atoi();
1125  TString variable = variables[num];
1126  TString pattern = TString::Format("{V%d}",num);
1127  i -= 2; // restore original position
1128  body.Replace(i, pattern.Length(),variable,variable.Length());
1129  i += variable.Length()-1; // update i to reflect change in body string
1130  }
1131  else if(body[i] == '[')
1132  {
1133  // update parameter counters in case of many functions (e.g. gaus(0)+gaus(3) )
1134  Int_t tmp = i;
1135  while(tmp < body.Length() && body[tmp] != ']')
1136  {
1137  tmp++;
1138  }
1139  Int_t num = TString(body(i+1,tmp - 1 - i)).Atoi();
1140  num += counter;
1141  TString replacement = TString::Format("%d",num);
1142 
1143  body.Replace(i+1,tmp - 1 - i,replacement,replacement.Length());
1144  i += replacement.Length() + 1;
1145  }
1146  }
1147  TString pattern;
1148  if(defaultCounter && defaultVariables)
1149  {
1150  pattern = TString::Format("%s%s",
1151  funName.Data(),
1152  (isNormalized ? "n" : ""));
1153  }
1154  if(!defaultCounter && defaultVariables)
1155  {
1156  pattern = TString::Format("%s%s(%d)",
1157  funName.Data(),
1158  (isNormalized ? "n" : ""),
1159  counter);
1160  }
1161  if(defaultCounter && !defaultVariables)
1162  {
1163  pattern = TString::Format("%s%s[%s]",
1164  funName.Data(),
1165  (isNormalized ? "n":""),
1166  varList.Data());
1167  }
1168  if(!defaultCounter && !defaultVariables)
1169  {
1170  pattern = TString::Format("%s%s[%s](%d)",
1171  funName.Data(),
1172  (isNormalized ? "n" : ""),
1173  varList.Data(),
1174  counter);
1175  }
1176  TString replacement = body;
1177 
1178  // set the number (only in case a function exists without anything else
1179  if (fNumber == 0 && formula.Length() <= (pattern.Length()-funPos) +1 ) { // leave 1 extra
1180  fNumber = functionsNumbers[funName] + 10*(dim-1);
1181  }
1182 
1183  //std::cout << " replace " << pattern << " with " << replacement << std::endl;
1184 
1185  formula.Replace(funPos,pattern.Length(),replacement,replacement.Length());
1186 
1187  funPos = formula.Index(funName);
1188  }
1189  //std::cout << " End loop of " << funName << " formula is now " << formula << std::endl;
1190  }
1191 
1192 }
1193 
1194 ////////////////////////////////////////////////////////////////////////////////
1195 /// Handling exponentiation
1196 /// Can handle multiple carets, eg.
1197 /// 2^3^4 will be treated like 2^(3^4)
1198 
1199 void TFormula::HandleExponentiation(TString &formula)
1200 {
1201  Int_t caretPos = formula.Last('^');
1202  while(caretPos != kNPOS && !IsAParameterName(formula,caretPos) )
1203  {
1204 
1205  TString right,left;
1206  Int_t temp = caretPos;
1207  temp--;
1208  // get the expression in ( ) which has the operator^ applied
1209  if(formula[temp] == ')')
1210  {
1211  Int_t depth = 1;
1212  temp--;
1213  while(depth != 0 && temp > 0)
1214  {
1215  if(formula[temp] == ')')
1216  depth++;
1217  if(formula[temp] == '(')
1218  depth--;
1219  temp--;
1220  }
1221  if (depth == 0) temp++;
1222  }
1223  // this in case of someting like sin(x+2)^2
1224  do {
1225  temp--; // go down one
1226  // handle scientific notation cases (1.e-2 ^ 3 )
1227  if (temp>=2 && IsScientificNotation(formula, temp-1) ) temp-=3;
1228  }
1229  while(temp >= 0 && !IsOperator(formula[temp]) && !IsBracket(formula[temp]) );
1230 
1231  assert(temp+1 >= 0);
1232  Int_t leftPos = temp+1;
1233  left = formula(leftPos, caretPos - leftPos);
1234  //std::cout << "left to replace is " << left << std::endl;
1235 
1236  // look now at the expression after the ^ operator
1237  temp = caretPos;
1238  temp++;
1239  if (temp >= formula.Length() ) {
1240  Error("HandleExponentiation","Invalid position of operator ^");
1241  return;
1242  }
1243  if(formula[temp] == '(')
1244  {
1245  Int_t depth = 1;
1246  temp++;
1247  while(depth != 0 && temp < formula.Length())
1248  {
1249  if(formula[temp] == ')')
1250  depth--;
1251  if(formula[temp] == '(')
1252  depth++;
1253  temp++;
1254  }
1255  temp--;
1256  }
1257  else {
1258  // handle case first character is operator - or + continue
1259  if (formula[temp] == '-' || formula[temp] == '+' ) temp++;
1260  // handle cases x^-2 or x^+2
1261  // need to handle also cases x^sin(x+y)
1262  Int_t depth = 0;
1263  // stop right expression if is an operator or if is a ")" from a zero depth
1264  while(temp < formula.Length() && ( (depth > 0) || !IsOperator(formula[temp]) ) )
1265  {
1266  temp++;
1267  // handle scientific notation cases (1.e-2 ^ 3 )
1268  if (temp>=2 && IsScientificNotation(formula, temp) ) temp+=2;
1269  // for internal parenthesis
1270  if (temp < formula.Length() && formula[temp] == '(') depth++;
1271  if (temp < formula.Length() && formula[temp] == ')') {
1272  if (depth > 0)
1273  depth--;
1274  else
1275  break; // case of end of a previously started expression e.g. sin(x^2)
1276  }
1277  }
1278  }
1279  right = formula(caretPos + 1, (temp - 1) - caretPos );
1280  //std::cout << "right to replace is " << right << std::endl;
1281 
1282  TString pattern = TString::Format("%s^%s",left.Data(),right.Data());
1283  TString replacement = TString::Format("pow(%s,%s)",left.Data(),right.Data());
1284 
1285  //std::cout << "pattern : " << pattern << std::endl;
1286  //std::cout << "replacement : " << replacement << std::endl;
1287  formula.Replace(leftPos,pattern.Length(),replacement,replacement.Length());
1288 
1289  caretPos = formula.Last('^');
1290  }
1291 }
1292 
1293 ////////////////////////////////////////////////////////////////////////////////
1294 /// Handle linear functions defined with the operator ++.
1295 
1296 void TFormula::HandleLinear(TString &formula)
1297 {
1298  // Handle Linear functions identified with "@" operator
1299  Int_t linPos = formula.Index("@");
1300  if (linPos == kNPOS ) return; // function is not linear
1301  Int_t nofLinParts = formula.CountChar((int)'@');
1302  assert(nofLinParts > 0);
1303  fLinearParts.reserve(nofLinParts + 1);
1304  Int_t Nlinear = 0;
1305  bool first = true;
1306  while(linPos != kNPOS && !IsAParameterName(formula,linPos))
1307  {
1308  SetBit(kLinear,1);
1309  // analyze left part only the first time
1310  Int_t temp = 0;
1311  TString left;
1312  if (first) {
1313  temp = linPos - 1;
1314  while(temp >= 0 && formula[temp] != '@')
1315  {
1316  temp--;
1317  }
1318  left = formula(temp+1,linPos - (temp +1));
1319  }
1320  temp = linPos + 1;
1321  while(temp < formula.Length() && formula[temp] != '@')
1322  {
1323  temp++;
1324  }
1325  TString right = formula(linPos+1,temp - (linPos+1));
1326 
1327  TString pattern = (first) ? TString::Format("%s@%s",left.Data(),right.Data()) : TString::Format("@%s",right.Data());
1328  TString replacement = (first) ? TString::Format("([%d]*(%s))+([%d]*(%s))",Nlinear,left.Data(),Nlinear+1,right.Data()) : TString::Format("+([%d]*(%s))",Nlinear,right.Data());
1329  Nlinear += (first) ? 2 : 1;
1330 
1331  formula.ReplaceAll(pattern,replacement);
1332  if (first) {
1333  TFormula *lin1 = new TFormula("__linear1",left,false);
1334  fLinearParts.push_back(lin1);
1335  }
1336  TFormula *lin2 = new TFormula("__linear2",right,false);
1337  fLinearParts.push_back(lin2);
1338 
1339  linPos = formula.Index("@");
1340  first = false;
1341  }
1342 }
1343 
1344 ////////////////////////////////////////////////////////////////////////////////
1345 /// Preprocessing of formula
1346 /// Replace all ** by ^, and removes spaces.
1347 /// Handle also parametrized functions like polN,gaus,expo,landau
1348 /// and exponentiation.
1349 /// Similar functionality should be added here.
1350 
1351 void TFormula::PreProcessFormula(TString &formula)
1352 {
1353  formula.ReplaceAll("**","^");
1354  formula.ReplaceAll("++","@"); // for linear functions
1355  formula.ReplaceAll(" ","");
1356  HandlePolN(formula);
1357  HandleParametrizedFunctions(formula);
1358  HandleExponentiation(formula);
1359  // "++" wil be dealt with Handle Linear
1360  HandleLinear(formula);
1361  // special case for "--" and "++"
1362  // ("++" needs to be written with whitespace that is removed before but then we re-add it again
1363  formula.ReplaceAll("--","- -");
1364  formula.ReplaceAll("++","+ +");
1365 }
1366 
1367 ////////////////////////////////////////////////////////////////////////////////
1368 /// prepare the formula to be executed
1369 /// normally is called with fFormula
1370 
1372 {
1373  fFuncs.clear();
1374  fReadyToExecute = false;
1375  ExtractFunctors(formula);
1376 
1377  // update the expression with the new formula
1378  fFormula = formula;
1379  // save formula to parse variable and parameters for Cling
1380  fClingInput = formula;
1381  // replace all { and }
1382  fFormula.ReplaceAll("{","");
1383  fFormula.ReplaceAll("}","");
1384 
1385  //std::cout << "functors are extracted formula is " << std::endl;
1386  //std::cout << fFormula << std::endl << std::endl;
1387 
1388  fFuncs.sort();
1389  fFuncs.unique();
1390 
1391  // use inputFormula for Cling
1393 
1394  // for pre-defined functions (need after processing)
1395  if (fNumber != 0) SetPredefinedParamNames();
1396 
1398 }
1399 
1400 ////////////////////////////////////////////////////////////////////////////////
1401 /// Extracts functors from formula, and put them in fFuncs.
1402 /// Simple grammar:
1403 /// - <function> := name(arg1,arg2...)
1404 /// - <variable> := name
1405 /// - <parameter> := [number]
1406 /// - <name> := String containing lower and upper letters, numbers, underscores
1407 /// - <number> := Integer number
1408 /// Operators are omitted.
1409 
1410 void TFormula::ExtractFunctors(TString &formula)
1411 {
1412  TString name = "";
1413  TString body = "";
1414  //printf("formula is : %s \n",formula.Data() );
1415  for(Int_t i = 0 ; i < formula.Length(); ++i )
1416  {
1417 
1418  //std::cout << "loop on character : " << i << " " << formula[i] << std::endl;
1419  // case of parameters
1420  if(formula[i] == '[')
1421  {
1422  Int_t tmp = i;
1423  i++;
1424  TString param = "";
1425  while(formula[i] != ']' && i < formula.Length())
1426  {
1427  param.Append(formula[i++]);
1428  }
1429  i++;
1430  //rename parameter name XX to pXX
1431  //std::cout << "examine parameters " << param << std::endl;
1432  int paramIndex = -1;
1433  if (param.IsDigit()) {
1434  paramIndex = param.Atoi();
1435  param.Insert(0, 'p'); // needed for the replacement
1436  if (paramIndex >= fNpar || fParams.find(param) == fParams.end()) {
1437  // add all parameters up to given index found
1438  for (int idx = 0; idx <= paramIndex; ++idx) {
1439  TString pname = TString::Format("p%d", idx);
1440  if (fParams.find(pname) == fParams.end())
1441  DoAddParameter(pname, 0, false);
1442  }
1443  }
1444  } else {
1445  // handle whitespace characters in parname
1446  param.ReplaceAll("\\s", " ");
1447  DoAddParameter(param, 0, false);
1448  }
1449  TString replacement = TString::Format("{[%s]}",param.Data());
1450  formula.Replace(tmp,i - tmp, replacement,replacement.Length());
1451  fFuncs.push_back(TFormulaFunction(param));
1452  //printf("found parameter %s \n",param.Data() );
1453  continue;
1454  }
1455  // case of strings
1456  if (formula[i] == '\"') {
1457  // look for next instance of "\"
1458  do {
1459  i++;
1460  } while(formula[i] != '\"');
1461  }
1462  // case of e or E for numbers in exponential notaton (e.g. 2.2e-3)
1463  if (IsScientificNotation(formula, i) )
1464  continue;
1465  // case of x for hexadecimal numbers
1466  if (IsHexadecimal(formula, i) ) {
1467  // find position of operator
1468  // do not check cases if character is not only a to f, but accept anything
1469  while ( !IsOperator(formula[i]) && i < formula.Length() ) {
1470  i++;
1471  }
1472  continue;
1473  }
1474 
1475 
1476  //std::cout << "investigating character : " << i << " " << formula[i] << " of formula " << formula << std::endl;
1477  // look for variable and function names. They start in C++ with alphanumeric characters
1478  if(isalpha(formula[i]) && !IsOperator(formula[i])) // not really needed to check if operator (if isalpha is not an operator)
1479  {
1480  //std::cout << "character : " << i << " " << formula[i] << " is not an operator and is alpha " << std::endl;
1481 
1482  while( IsFunctionNameChar(formula[i]) && i < formula.Length())
1483  {
1484  // need special case for separting operator ":" from scope operator "::"
1485  if (formula[i] == ':' && ( (i+1) < formula.Length() ) ) {
1486  if ( formula[i+1] == ':' ) {
1487  // case of :: (scopeOperator)
1488  name.Append("::");
1489  i+=2;
1490  continue;
1491  }
1492  else
1493  break;
1494  }
1495 
1496  name.Append(formula[i++]);
1497  }
1498  //printf(" build a name %s \n",name.Data() );
1499  if(formula[i] == '(')
1500  {
1501  i++;
1502  if(formula[i] == ')')
1503  {
1504  fFuncs.push_back(TFormulaFunction(name,body,0));
1505  name = body = "";
1506  continue;
1507  }
1508  Int_t depth = 1;
1509  Int_t args = 1; // we will miss first argument
1510  while(depth != 0 && i < formula.Length())
1511  {
1512  switch(formula[i])
1513  {
1514  case '(': depth++; break;
1515  case ')': depth--; break;
1516  case ',': if(depth == 1) args++; break;
1517  }
1518  if(depth != 0) // we don't want last ')' inside body
1519  {
1520  body.Append(formula[i++]);
1521  }
1522  }
1523  Int_t originalBodyLen = body.Length();
1524  ExtractFunctors(body);
1525  formula.Replace(i-originalBodyLen,originalBodyLen,body,body.Length());
1526  i += body.Length() - originalBodyLen;
1527  fFuncs.push_back(TFormulaFunction(name,body,args));
1528  }
1529  else
1530  {
1531 
1532  //std::cout << "check if character : " << i << " " << formula[i] << " from name " << name << " is a function " << std::endl;
1533 
1534  // check if function is provided by gROOT
1535  TObject *obj = 0;
1536  {
1538  obj = gROOT->GetListOfFunctions()->FindObject(name);
1539  }
1540  TFormula * f = dynamic_cast<TFormula*> (obj);
1541  if (!f) {
1542  // maybe object is a TF1
1543  TF1 * f1 = dynamic_cast<TF1*> (obj);
1544  if (f1) f = f1->GetFormula();
1545  }
1546  if (f) {
1547  TString replacementFormula = f->GetExpFormula();
1548  // analyze expression string
1549  //std::cout << "formula to replace for " << f->GetName() << " is " << replacementFormula << std::endl;
1550  PreProcessFormula(replacementFormula);
1551  // we need to define different parameters if we use the unnamed default parameters ([0])
1552  // I need to replace all the terms in the functor for backward compatibility of the case
1553  // f1("[0]*x") f2("[0]*x") f1+f2 - it is weird but it is better to support
1554  //std::cout << "current number of parameter is " << fNpar << std::endl;
1555  int nparOffset = 0;
1556  //if (fParams.find("0") != fParams.end() ) {
1557  // do in any case if parameters are existing
1558  std::vector<TString> newNames;
1559  if (fNpar > 0) {
1560  nparOffset = fNpar;
1561  newNames.resize(f->GetNpar() );
1562  // start from higher number to avoid overlap
1563  for (int jpar = f->GetNpar()-1; jpar >= 0; --jpar ) {
1564  // parameters name have a "p" added in front
1565  TString pj = TString(f->GetParName(jpar));
1566  if ( pj[0] == 'p' && TString(pj(1,pj.Length())).IsDigit() ) {
1567  TString oldName = TString::Format("[%s]",f->GetParName(jpar));
1568  TString newName = TString::Format("[p%d]",nparOffset+jpar);
1569  //std::cout << "replace - parameter " << f->GetParName(jpar) << " with " << newName << std::endl;
1570  replacementFormula.ReplaceAll(oldName,newName);
1571  newNames[jpar] = newName;
1572  }
1573  else
1574  newNames[jpar] = f->GetParName(jpar);
1575  }
1576  //std::cout << "after replacing params " << replacementFormula << std::endl;
1577  }
1578  ExtractFunctors(replacementFormula);
1579  //std::cout << "after re-extracting functors " << replacementFormula << std::endl;
1580 
1581  // set parameter value from replacement formula
1582  for (int jpar = 0; jpar < f->GetNpar(); ++jpar) {
1583  if (nparOffset> 0) {
1584  // parameter have an offset- so take this into accound
1585  assert((int) newNames.size() == f->GetNpar() );
1586  SetParameter(newNames[jpar], f->GetParameter(jpar) );
1587  }
1588  else
1589  // names are the same between current formula and replaced one
1590  SetParameter(f->GetParName(jpar), f->GetParameter(jpar) );
1591  }
1592  // need to add parenthesis at begin end end of replacementFormula
1593  replacementFormula.Insert(0,'(');
1594  replacementFormula.Insert(replacementFormula.Length(),')');
1595  formula.Replace(i-name.Length(),name.Length(), replacementFormula, replacementFormula.Length());
1596  // move forward the index i of the main loop
1597  i += replacementFormula.Length()-name.Length();
1598 
1599  // we have extracted all the functor for "fname"
1600  //std::cout << " i = " << i << " f[i] = " << formula[i] << " - " << formula << std::endl;
1601  name = "";
1602 
1603  continue;
1604  }
1605 
1606  // add now functor in
1607  TString replacement = TString::Format("{%s}",name.Data());
1608  formula.Replace(i-name.Length(),name.Length(),replacement,replacement.Length());
1609  i += 2;
1610  fFuncs.push_back(TFormulaFunction(name));
1611  }
1612  }
1613  name = body = "";
1614 
1615  }
1616 }
1617 
1618 ////////////////////////////////////////////////////////////////////////////////
1619 /// Iterates through funtors in fFuncs and performs the appropriate action.
1620 /// If functor has 0 arguments (has only name) can be:
1621 /// - variable
1622 /// * will be replaced with x[num], where x is an array containing value of this variable under num.
1623 /// - pre-defined formula
1624 /// * will be replaced with formulas body
1625 /// - constant
1626 /// * will be replaced with constant value
1627 /// - parameter
1628 /// * will be replaced with p[num], where p is an array containing value of this parameter under num.
1629 /// If has arguments it can be :
1630 /// - function shortcut, eg. sin
1631 /// * will be replaced with fullname of function, eg. sin -> TMath::Sin
1632 /// - function from cling environment, eg. TMath::BreitWigner(x,y,z)
1633 /// * first check if function exists, and has same number of arguments, then accept it and set as found.
1634 /// If all functors after iteration are matched with corresponding action,
1635 /// it inputs C++ code of formula into cling, and sets flag that formula is ready to evaluate.
1636 
1637 void TFormula::ProcessFormula(TString &formula)
1638 {
1639  // std::cout << "Begin: formula is " << formula << " list of functors " << fFuncs.size() << std::endl;
1640 
1641  for(list<TFormulaFunction>::iterator funcsIt = fFuncs.begin(); funcsIt != fFuncs.end(); ++funcsIt)
1642  {
1643  TFormulaFunction & fun = *funcsIt;
1644 
1645  //std::cout << "fun is " << fun.GetName() << std::endl;
1646 
1647  if(fun.fFound)
1648  continue;
1649  if(fun.IsFuncCall())
1650  {
1651  map<TString,TString>::iterator it = fFunctionsShortcuts.find(fun.GetName());
1652  if(it != fFunctionsShortcuts.end())
1653  {
1654  TString shortcut = it->first;
1655  TString full = it->second;
1656  //std::cout << " functor " << fun.GetName() << " found - replace " << shortcut << " with " << full << " in " << formula << std::endl;
1657  // replace all functors
1658  Ssiz_t index = formula.Index(shortcut,0);
1659  while ( index != kNPOS) {
1660  // check that function is not in a namespace and is not in other characters
1661  //std::cout << "analyzing " << shortcut << " in " << formula << std::endl;
1662  Ssiz_t i2 = index + shortcut.Length();
1663  if ( (index > 0) && (isalpha( formula[index-1] ) || formula[index-1] == ':' )) {
1664  index = formula.Index(shortcut,i2);
1665  continue;
1666  }
1667  if (i2 < formula.Length() && formula[i2] != '(') {
1668  index = formula.Index(shortcut,i2);
1669  continue;
1670  }
1671  // now replace the string
1672  formula.Replace(index, shortcut.Length(), full);
1673  Ssiz_t inext = index + full.Length();
1674  index = formula.Index(shortcut,inext);
1675  fun.fFound = true;
1676  }
1677  }
1678  if(fun.fName.Contains("::")) // add support for nested namespaces
1679  {
1680  // look for last occurence of "::"
1681  std::string name(fun.fName);
1682  size_t index = name.rfind("::");
1683  assert(index != std::string::npos);
1684  TString className = fun.fName(0,fun.fName(0,index).Length());
1685  TString functionName = fun.fName(index + 2, fun.fName.Length());
1686 
1687  Bool_t silent = true;
1688  TClass *tclass = TClass::GetClass(className,silent);
1689  // std::cout << "looking for class " << className << std::endl;
1690  const TList *methodList = tclass->GetListOfAllPublicMethods();
1691  TIter next(methodList);
1692  TMethod *p;
1693  while ((p = (TMethod*) next()))
1694  {
1695  if (strcmp(p->GetName(),functionName.Data()) == 0 &&
1696  (fun.GetNargs() <= p->GetNargs() && fun.GetNargs() >= p->GetNargs() - p->GetNargsOpt() ) )
1697  {
1698  fun.fFound = true;
1699  break;
1700  }
1701  }
1702  }
1703  if(!fun.fFound)
1704  {
1705  // try to look into all the global functions in gROOT
1706  TFunction* f;
1707  {
1709  f = (TFunction*) gROOT->GetListOfGlobalFunctions(true)->FindObject(fun.fName);
1710  }
1711  // if found a function with matching arguments
1712  if (f && fun.GetNargs() <= f->GetNargs() && fun.GetNargs() >= f->GetNargs() - f->GetNargsOpt() )
1713  {
1714  fun.fFound = true;
1715  }
1716  }
1717 
1718  if(!fun.fFound)
1719  {
1720  // ignore not found functions
1721  if (gDebug)
1722  Info("TFormula","Could not find %s function with %d argument(s)",fun.GetName(),fun.GetNargs());
1723  fun.fFound = false;
1724  }
1725  }
1726  else
1727  {
1728  TFormula* old = 0;
1729  {
1731  old = (TFormula*)gROOT->GetListOfFunctions()->FindObject(gNamePrefix + fun.fName);
1732  }
1733  if(old)
1734  {
1735  // we should not go here (this analysis is done before in ExtractFunctors)
1736  assert(false);
1737  fun.fFound = true;
1738  TString pattern = TString::Format("{%s}",fun.GetName());
1739  TString replacement = old->GetExpFormula();
1740  PreProcessFormula(replacement);
1741  ExtractFunctors(replacement);
1742  formula.ReplaceAll(pattern,replacement);
1743  continue;
1744  }
1745  // looking for default variables defined in fVars
1746 
1747  map<TString,TFormulaVariable>::iterator varsIt = fVars.find(fun.GetName());
1748  if(varsIt!= fVars.end())
1749  {
1750 
1751  TString name = (*varsIt).second.GetName();
1752  Double_t value = (*varsIt).second.fValue;
1753 
1754 
1755  AddVariable(name,value); // this set the cling variable
1756  if(!fVars[name].fFound)
1757  {
1758 
1759 
1760  fVars[name].fFound = true;
1761  int varDim = (*varsIt).second.fArrayPos; // variable dimenions (0 for x, 1 for y, 2, for z)
1762  if (varDim >= fNdim) {
1763  fNdim = varDim+1;
1764 
1765  // we need to be sure that all other variables are added with position less
1766  for ( auto &v : fVars) {
1767  if (v.second.fArrayPos < varDim && !v.second.fFound ) {
1768  AddVariable(v.first, v.second.fValue);
1769  v.second.fFound = true;
1770  }
1771  }
1772  }
1773  }
1774  // remove the "{.. }" added around the variable
1775  TString pattern = TString::Format("{%s}",name.Data());
1776  TString replacement = TString::Format("x[%d]",(*varsIt).second.fArrayPos);
1777  formula.ReplaceAll(pattern,replacement);
1778 
1779  //std::cout << "Found an observable for " << fun.GetName() << std::endl;
1780 
1781  fun.fFound = true;
1782  continue;
1783  }
1784  // check for observables defined as x[0],x[1],....
1785  // maybe could use a regular expression here
1786  // only in case match with defined variables is not successfull
1787  TString funname = fun.GetName();
1788  if (funname.Contains("x[") && funname.Contains("]") ) {
1789  TString sdigit = funname(2,funname.Index("]") );
1790  int digit = sdigit.Atoi();
1791  if (digit >= fNdim) {
1792  fNdim = digit+1;
1793  // we need to add the variables in fVars all of them before x[n]
1794  for (int j = 0; j < fNdim; ++j) {
1795  TString vname = TString::Format("x[%d]",j);
1796  if (fVars.find(vname) == fVars.end() ) {
1797  fVars[vname] = TFormulaVariable(vname,0,j);
1798  fVars[vname].fFound = true;
1799  AddVariable(vname,0.);
1800  }
1801  }
1802  }
1803  //std::cout << "Found matching observable for " << funname << std::endl;
1804  fun.fFound = true;
1805  // remove the "{.. }" added around the variable
1806  TString pattern = TString::Format("{%s}",funname.Data());
1807  formula.ReplaceAll(pattern,funname);
1808  continue;
1809  }
1810  //}
1811 
1812  auto paramsIt = fParams.find(fun.GetName());
1813  if(paramsIt != fParams.end())
1814  {
1815  //TString name = (*paramsIt).second.GetName();
1816  TString pattern = TString::Format("{[%s]}",fun.GetName());
1817  //std::cout << "pattern is " << pattern << std::endl;
1818  if(formula.Index(pattern) != kNPOS)
1819  {
1820  //TString replacement = TString::Format("p[%d]",(*paramsIt).second.fArrayPos);
1821  TString replacement = TString::Format("p[%d]",(*paramsIt).second);
1822  //std::cout << "replace pattern " << pattern << " with " << replacement << std::endl;
1823  formula.ReplaceAll(pattern,replacement);
1824 
1825  }
1826  fun.fFound = true;
1827  continue;
1828  }
1829  else {
1830  //std::cout << "functor " << fun.GetName() << " is not a parameter " << std::endl;
1831  }
1832 
1833  // looking for constants (needs to be done after looking at the parameters)
1834  map<TString,Double_t>::iterator constIt = fConsts.find(fun.GetName());
1835  if(constIt != fConsts.end())
1836  {
1837  TString pattern = TString::Format("{%s}",fun.GetName());
1838  TString value = TString::Format("%lf",(*constIt).second);
1839  formula.ReplaceAll(pattern,value);
1840  fun.fFound = true;
1841  //std::cout << "constant with name " << fun.GetName() << " is found " << std::endl;
1842  continue;
1843  }
1844 
1845 
1846  fun.fFound = false;
1847  }
1848  }
1849  //std::cout << "End: formula is " << formula << std::endl;
1850 
1851  // ignore case of functors have been matched - try to pass it to Cling
1852  if(!fReadyToExecute)
1853  {
1854  fReadyToExecute = true;
1855  Bool_t hasVariables = (fNdim > 0);
1856  Bool_t hasParameters = (fNpar > 0);
1857  if(!hasParameters)
1858  {
1859  fAllParametersSetted = true;
1860  }
1861  // assume a function without variables is always 1-dimensional
1862  if (hasParameters && ! hasVariables) {
1863  fNdim = 1;
1864  AddVariable("x",0);
1865  hasVariables = true;
1866  }
1867  Bool_t hasBoth = hasVariables && hasParameters;
1868  Bool_t inputIntoCling = (formula.Length() > 0);
1869  if (inputIntoCling) {
1870 
1871  // save copy of inputFormula in a std::strig for the unordered map
1872  // and also formula is same as FClingInput typically and it will be modified
1873  std::string inputFormula = std::string(formula);
1874 
1875 
1876  // valid input formula - try to put into Cling
1877  TString argumentsPrototype =
1878  TString::Format("%s%s%s",(hasVariables ? "Double_t *x" : ""), (hasBoth ? "," : ""),
1879  (hasParameters ? "Double_t *p" : ""));
1880 
1881 
1882  // set the name for Cling using the hash_function
1884 
1885  // check if formula exist already in the map
1887 
1888  auto funcit = gClingFunctions.find(inputFormula);
1889 
1890  if (funcit != gClingFunctions.end() ) {
1892  fClingInitialized = true;
1893  inputIntoCling = false;
1894  }
1895 
1896  // set the cling name using hash of the static formulae map
1897  auto hasher = gClingFunctions.hash_function();
1898  fClingName = TString::Format("%s__id%zu",gNamePrefix.Data(), hasher(inputFormula) );
1899 
1900  fClingInput = TString::Format("Double_t %s(%s){ return %s ; }", fClingName.Data(),argumentsPrototype.Data(),inputFormula.c_str());
1901 
1902  // this is not needed (maybe can be re-added in case of recompilation of identical expressions
1903  // // check in case of a change if need to re-initialize
1904  // if (fClingInitialized) {
1905  // if (oldClingInput == fClingInput)
1906  // inputIntoCling = false;
1907  // else
1908  // fClingInitialized = false;
1909  // }
1910 
1911 
1912  if(inputIntoCling) {
1914  if (fClingInitialized) {
1915  // if Cling has been succesfully initialized
1916  // dave function ptr in the static map
1918  gClingFunctions.insert ( std::make_pair ( inputFormula, (void*) fFuncPtr) );
1919  }
1920 
1921  }
1922  else {
1923  fAllParametersSetted = true;
1924  fClingInitialized = true;
1925  }
1926  }
1927  }
1928 
1929 
1930  // IN case of a Cling Error check components wich are not found in Cling
1931  // check that all formula components arematched otherwise emit an error
1932  if (!fClingInitialized) {
1933  Bool_t allFunctorsMatched = true;
1934  for(list<TFormulaFunction>::iterator it = fFuncs.begin(); it != fFuncs.end(); it++)
1935  {
1936  if(!it->fFound)
1937  {
1938  allFunctorsMatched = false;
1939  if (it->GetNargs() == 0)
1940  Error("ProcessFormula","\"%s\" has not been matched in the formula expression",it->GetName() );
1941  else
1942  Error("ProcessFormula","Could not find %s function with %d argument(s)",it->GetName(),it->GetNargs());
1943  }
1944  }
1945  if (!allFunctorsMatched) {
1946  Error("ProcessFormula","Formula \"%s\" is invalid !", GetExpFormula().Data() );
1947  fReadyToExecute = false;
1948  }
1949  }
1950 
1951  // clean up un-used default variables in case formula is valid
1953  auto itvar = fVars.begin();
1954  do
1955  {
1956  if ( ! itvar->second.fFound ) {
1957  //std::cout << "Erase variable " << itvar->first << std::endl;
1958  itvar = fVars.erase(itvar);
1959  }
1960  else
1961  itvar++;
1962  }
1963  while( itvar != fVars.end() );
1964  }
1965 
1966 }
1967 
1968 ////////////////////////////////////////////////////////////////////////////////
1969 /// Set parameter names only in case of pre-defined functions.
1970 
1972 
1973  if (fNumber == 0) return;
1974 
1975  if (fNumber == 100) { // Gaussian
1976  SetParName(0,"Constant");
1977  SetParName(1,"Mean");
1978  SetParName(2,"Sigma");
1979  return;
1980  }
1981  if (fNumber == 110) {
1982  SetParName(0,"Constant");
1983  SetParName(1,"MeanX");
1984  SetParName(2,"SigmaX");
1985  SetParName(3,"MeanY");
1986  SetParName(4,"SigmaY");
1987  return;
1988  }
1989  if (fNumber == 112) { // bigaus
1990  SetParName(0,"Constant");
1991  SetParName(1,"MeanX");
1992  SetParName(2,"SigmaX");
1993  SetParName(3,"MeanY");
1994  SetParName(4,"SigmaY");
1995  SetParName(5,"Rho");
1996  return;
1997  }
1998  if (fNumber == 200) { // exponential
1999  SetParName(0,"Constant");
2000  SetParName(1,"Slope");
2001  return;
2002  }
2003  if (fNumber == 400) { // landau
2004  SetParName(0,"Constant");
2005  SetParName(1,"MPV");
2006  SetParName(2,"Sigma");
2007  return;
2008  }
2009  if (fNumber == 500) { // crystal-ball
2010  SetParName(0,"Constant");
2011  SetParName(1,"Mean");
2012  SetParName(2,"Sigma");
2013  SetParName(3,"Alpha");
2014  SetParName(4,"N");
2015  return;
2016  }
2017  if (fNumber == 600) { // breit-wigner
2018  SetParName(0,"Constant");
2019  SetParName(1,"Mean");
2020  SetParName(2,"Gamma");
2021  return;
2022  }
2023  // if formula is a polynomial (or chebyshev), set parameter names
2024  // not needed anymore (p0 is assigned by default)
2025  // if (fNumber == (300+fNpar-1) ) {
2026  // for (int i = 0; i < fNpar; i++) SetParName(i,TString::Format("p%d",i));
2027  // return;
2028  // }
2029 
2030  // // general case if parameters are digits (XX) change to pXX
2031  // auto paramMap = fParams; // need to copy the map because SetParName is going to modify it
2032  // for ( auto & p : paramMap) {
2033  // if (p.first.IsDigit() )
2034  // SetParName(p.second,TString::Format("p%s",p.first.Data()));
2035  // }
2036 
2037  return;
2038 }
2039 
2040 ////////////////////////////////////////////////////////////////////////////////
2041 /// Return linear part.
2042 
2044 {
2045  if (!fLinearParts.empty()) {
2046  int n = fLinearParts.size();
2047  if (i < 0 || i >= n ) {
2048  Error("GetLinearPart","Formula %s has only %d linear parts - requested %d",GetName(),n,i);
2049  return nullptr;
2050  }
2051  return fLinearParts[i];
2052  }
2053  return nullptr;
2054 }
2055 
2056 ////////////////////////////////////////////////////////////////////////////////
2057 /// Adds variable to known variables, and reprocess formula.
2058 
2059 void TFormula::AddVariable(const TString &name, double value)
2060 {
2061  if(fVars.find(name) != fVars.end() )
2062  {
2063  TFormulaVariable & var = fVars[name];
2064  var.fValue = value;
2065 
2066  // If the position is not defined in the Cling vectors, make space for it
2067  // but normally is variable is defined in fVars a slot should be also present in fClingVariables
2068  if(var.fArrayPos < 0)
2069  {
2070  var.fArrayPos = fVars.size();
2071  }
2072  if(var.fArrayPos >= (int)fClingVariables.size())
2073  {
2074  fClingVariables.resize(var.fArrayPos+1);
2075  }
2076  fClingVariables[var.fArrayPos] = value;
2077  }
2078  else
2079  {
2080  TFormulaVariable var(name,value,fVars.size());
2081  fVars[name] = var;
2082  fClingVariables.push_back(value);
2083  if (!fFormula.IsNull() ) {
2084  //printf("process formula again - %s \n",fClingInput.Data() );
2086  }
2087  }
2088 
2089 }
2090 
2091 ////////////////////////////////////////////////////////////////////////////////
2092 /// Adds multiple variables.
2093 /// First argument is an array of pairs<TString,Double>, where
2094 /// first argument is name of variable,
2095 /// second argument represents value.
2096 /// size - number of variables passed in first argument
2097 
2098 void TFormula::AddVariables(const TString *vars, const Int_t size)
2099 {
2100  Bool_t anyNewVar = false;
2101  for(Int_t i = 0 ; i < size; ++i)
2102  {
2103 
2104  const TString & vname = vars[i];
2105 
2106  TFormulaVariable &var = fVars[vname];
2107  if(var.fArrayPos < 0)
2108  {
2109 
2110  var.fName = vname;
2111  var.fArrayPos = fVars.size();
2112  anyNewVar = true;
2113  var.fValue = 0;
2114  if(var.fArrayPos >= (int)fClingVariables.capacity())
2115  {
2116  Int_t multiplier = 2;
2117  if(fFuncs.size() > 100)
2118  {
2119  multiplier = TMath::Floor(TMath::Log10(fFuncs.size()) * 10);
2120  }
2121  fClingVariables.reserve(multiplier * fClingVariables.capacity());
2122  }
2123  fClingVariables.push_back(0.0);
2124  }
2125  // else
2126  // {
2127  // var.fValue = v.second;
2128  // fClingVariables[var.fArrayPos] = v.second;
2129  // }
2130  }
2131  if(anyNewVar && !fFormula.IsNull())
2132  {
2134  }
2135 
2136 }
2137 
2138 ////////////////////////////////////////////////////////////////////////////////
2139 /// Set the name of the formula. We need to allow the list of function to
2140 /// properly handle the hashes.
2141 
2142 void TFormula::SetName(const char* name)
2143 {
2144  if (IsReservedName(name)) {
2145  Error("SetName","The name \'%s\' is reserved as a TFormula variable name.\n"
2146  "\tThis function will not be renamed.",name);
2147  } else {
2148  // Here we need to remove and re-add to keep the hashes consistent with
2149  // the underlying names.
2150  auto listOfFunctions = gROOT->GetListOfFunctions();
2151  TObject* thisAsFunctionInList = nullptr;
2153  if (listOfFunctions){
2154  thisAsFunctionInList = listOfFunctions->FindObject(this);
2155  if (thisAsFunctionInList) listOfFunctions->Remove(thisAsFunctionInList);
2156  }
2157  TNamed::SetName(name);
2158  if (thisAsFunctionInList) listOfFunctions->Add(thisAsFunctionInList);
2159  }
2160 }
2161 
2162 ////////////////////////////////////////////////////////////////////////////////
2163 ///
2164 /// Sets multiple variables.
2165 /// First argument is an array of pairs<TString,Double>, where
2166 /// first argument is name of variable,
2167 /// second argument represents value.
2168 /// size - number of variables passed in first argument
2169 
2170 void TFormula::SetVariables(const pair<TString,Double_t> *vars, const Int_t size)
2171 {
2172  for(Int_t i = 0; i < size; ++i)
2173  {
2174  pair<TString,Double_t> v = vars[i];
2175  if(fVars.find(v.first) != fVars.end())
2176  {
2177  fVars[v.first].fValue = v.second;
2178  fClingVariables[fVars[v.first].fArrayPos] = v.second;
2179  }
2180  else
2181  {
2182  Error("SetVariables","Variable %s is not defined.",v.first.Data());
2183  }
2184  }
2185 }
2186 
2187 ////////////////////////////////////////////////////////////////////////////////
2188 /// Returns variable value.
2189 
2191 {
2192  TString sname(name);
2193  if(fVars.find(sname) == fVars.end())
2194  {
2195  Error("GetVariable","Variable %s is not defined.",sname.Data());
2196  return -1;
2197  }
2198  return fVars.find(sname)->second.fValue;
2199 }
2200 
2201 ////////////////////////////////////////////////////////////////////////////////
2202 /// Returns variable number (positon in array) given its name.
2203 
2205 {
2206  TString sname(name);
2207  if(fVars.find(sname) == fVars.end())
2208  {
2209  Error("GetVarNumber","Variable %s is not defined.",sname.Data());
2210  return -1;
2211  }
2212  return fVars.find(sname)->second.fArrayPos;
2213 }
2214 
2215 ////////////////////////////////////////////////////////////////////////////////
2216 /// Returns variable name given its position in the array.
2217 
2218 TString TFormula::GetVarName(Int_t ivar) const
2219 {
2220  if (ivar < 0 || ivar >= fNdim) return "";
2221 
2222  // need to loop on the map to find corresponding variable
2223  for ( auto & v : fVars) {
2224  if (v.second.fArrayPos == ivar) return v.first;
2225  }
2226  Error("GetVarName","Variable with index %d not found !!",ivar);
2227  //return TString::Format("x%d",ivar);
2228  return TString();
2229 }
2230 
2231 ////////////////////////////////////////////////////////////////////////////////
2232 /// Sets variable value.
2233 
2234 void TFormula::SetVariable(const TString &name, Double_t value)
2235 {
2236  if(fVars.find(name) == fVars.end())
2237  {
2238  Error("SetVariable","Variable %s is not defined.",name.Data());
2239  return;
2240  }
2241  fVars[name].fValue = value;
2242  fClingVariables[fVars[name].fArrayPos] = value;
2243 }
2244 
2245 ////////////////////////////////////////////////////////////////////////////////
2246 /// Adds parameter to known parameters.
2247 /// User should use SetParameter, because parameters are added during initialization part,
2248 /// and after that adding new will be pointless.
2249 
2250 void TFormula::DoAddParameter(const TString &name, Double_t value, Bool_t processFormula)
2251 {
2252  //std::cout << "adding parameter " << name << std::endl;
2253 
2254  // if parameter is already defined in fParams - just set the new value
2255  if(fParams.find(name) != fParams.end() )
2256  {
2257  int ipos = fParams[name];
2258  //TFormulaVariable & par = fParams[name];
2259  //par.fValue = value;
2260  if (ipos < 0)
2261  {
2262  ipos = fParams.size();
2263  fParams[name] = ipos;
2264  }
2265 //
2266  if(ipos >= (int)fClingParameters.size())
2267  {
2268  if(ipos >= (int)fClingParameters.capacity())
2269  fClingParameters.reserve( TMath::Max(int(fParams.size()), ipos+1));
2270  fClingParameters.insert(fClingParameters.end(),ipos+1-fClingParameters.size(),0.0);
2271  }
2272  fClingParameters[ipos] = value;
2273  }
2274  else
2275  {
2276  // new parameter defined
2277  fNpar++;
2278  //TFormulaVariable(name,value,fParams.size());
2279  int pos = fParams.size();
2280  //fParams.insert(std::make_pair<TString,TFormulaVariable>(name,TFormulaVariable(name,value,pos)));
2281  auto ret = fParams.insert(std::make_pair(name,pos));
2282  // map returns a std::pair<iterator, bool>
2283  // use map the order for defult position of parameters in the vector
2284  // (i.e use the alphabetic order)
2285  if (ret.second) {
2286  // a new element is inserted
2287  if (ret.first == fParams.begin() )
2288  pos = 0;
2289  else {
2290  auto previous = (ret.first);
2291  --previous;
2292  pos = previous->second + 1;
2293  }
2294 
2295 
2296  if (pos < (int)fClingParameters.size() )
2297  fClingParameters.insert(fClingParameters.begin()+pos,value);
2298  else {
2299  // this should not happen
2300  if (pos > (int)fClingParameters.size() )
2301  Warning("inserting parameter %s at pos %d when vector size is %d \n",name.Data(),pos,(int)fClingParameters.size() );
2302 
2303  if(pos >= (int)fClingParameters.capacity())
2304  fClingParameters.reserve( TMath::Max(int(fParams.size()), pos+1));
2305  fClingParameters.insert(fClingParameters.end(),pos+1-fClingParameters.size(),0.0);
2306  fClingParameters[pos] = value;
2307  }
2308 
2309  // need to adjust all other positions
2310  for ( auto it = ret.first; it != fParams.end(); ++it ) {
2311  it->second = pos;
2312  pos++;
2313  }
2314  // for (auto & p : fParams)
2315  // std::cout << "Parameter " << p.first << " position " << p.second << std::endl;
2316  // printf("inserted parameters size params %d size cling %d \n",fParams.size(), fClingParameters.size() );
2317  }
2318  if (processFormula) {
2319  // replace first in input parameter name with [name]
2320  fClingInput.ReplaceAll(name,TString::Format("[%s]",name.Data() ) );
2322  }
2323  }
2324 
2325 }
2326 
2327 ////////////////////////////////////////////////////////////////////////////////
2328 /// Return parameter index given a name (return -1 for not existing parameters)
2329 /// non need to print an error
2330 
2331 Int_t TFormula::GetParNumber(const char * name) const {
2332  auto it = fParams.find(name);
2333  if(it == fParams.end())
2334  {
2335  return -1;
2336  }
2337  return it->second;
2338 
2339 }
2340 
2341 ////////////////////////////////////////////////////////////////////////////////
2342 /// Returns parameter value given by string.
2343 
2345 {
2346  int i = GetParNumber(name);
2347  if (i == -1) {
2348  Error("GetParameter","Parameter %s is not defined.",name);
2349  return TMath::QuietNaN();
2350  }
2351 
2352  return GetParameter( GetParNumber(name) );
2353 }
2354 
2355 ////////////////////////////////////////////////////////////////////////////////
2356 /// Return parameter value given by integer.
2357 
2359 {
2360  //TString name = TString::Format("%d",param);
2361  if(param >=0 && param < (int) fClingParameters.size())
2362  return fClingParameters[param];
2363  Error("GetParameter","wrong index used - use GetParameter(name)");
2364  return TMath::QuietNaN();
2365 }
2366 
2367 ////////////////////////////////////////////////////////////////////////////////
2368 /// Return parameter name given by integer.
2369 
2370 const char * TFormula::GetParName(Int_t ipar) const
2371 {
2372  if (ipar < 0 || ipar >= fNpar) return "";
2373 
2374  // need to loop on the map to find corresponding parameter
2375  for ( auto & p : fParams) {
2376  if (p.second == ipar) return p.first.Data();
2377  }
2378  Error("GetParName","Parameter with index %d not found !!",ipar);
2379  //return TString::Format("p%d",ipar);
2380  return TString();
2381 }
2382 
2383 ////////////////////////////////////////////////////////////////////////////////
2385 {
2386  if(!fClingParameters.empty())
2387  return const_cast<Double_t*>(&fClingParameters[0]);
2388  return 0;
2389 }
2390 
2392 {
2393  for(Int_t i = 0; i < fNpar; ++i)
2394  {
2395  if (Int_t(fClingParameters.size()) > i)
2396  params[i] = fClingParameters[i];
2397  else
2398  params[i] = -1;
2399  }
2400 }
2401 
2402 ////////////////////////////////////////////////////////////////////////////////
2403 /// Sets parameter value.
2404 
2405 void TFormula::SetParameter(const char *name, Double_t value)
2406 {
2407  SetParameter( GetParNumber(name), value);
2408 
2409  // do we need this ???
2410 #ifdef OLDPARAMS
2411  if(fParams.find(name) == fParams.end())
2412  {
2413  Error("SetParameter","Parameter %s is not defined.",name.Data());
2414  return;
2415  }
2416  fParams[name].fValue = value;
2417  fParams[name].fFound = true;
2418  fClingParameters[fParams[name].fArrayPos] = value;
2419  fAllParametersSetted = true;
2420  for(map<TString,TFormulaVariable>::iterator it = fParams.begin(); it != fParams.end(); it++)
2421  {
2422  if(!it->second.fFound)
2423  {
2424  fAllParametersSetted = false;
2425  break;
2426  }
2427  }
2428 #endif
2429 }
2430 
2431 #ifdef OLDPARAMS
2432 
2433 ////////////////////////////////////////////////////////////////////////////////
2434 /// Set multiple parameters.
2435 /// First argument is an array of pairs<TString,Double>, where
2436 /// first argument is name of parameter,
2437 /// second argument represents value.
2438 /// size - number of params passed in first argument
2439 
2440 void TFormula::SetParameters(const pair<TString,Double_t> *params,const Int_t size)
2441 {
2442  for(Int_t i = 0 ; i < size ; ++i)
2443  {
2444  pair<TString,Double_t> p = params[i];
2445  if(fParams.find(p.first) == fParams.end())
2446  {
2447  Error("SetParameters","Parameter %s is not defined",p.first.Data());
2448  continue;
2449  }
2450  fParams[p.first].fValue = p.second;
2451  fParams[p.first].fFound = true;
2452  fClingParameters[fParams[p.first].fArrayPos] = p.second;
2453  }
2454  fAllParametersSetted = true;
2455  for(map<TString,TFormulaVariable>::iterator it = fParams.begin(); it != fParams.end(); it++)
2456  {
2457  if(!it->second.fFound)
2458  {
2459  fAllParametersSetted = false;
2460  break;
2461  }
2462  }
2463 }
2464 #endif
2465 
2466 ////////////////////////////////////////////////////////////////////////////////
2467 void TFormula::DoSetParameters(const Double_t *params, Int_t size)
2468 {
2469  if(!params || size < 0 || size > fNpar) return;
2470  // reset vector of cling parameters
2471  if (size != (int) fClingParameters.size() ) {
2472  Warning("SetParameters","size is not same of cling parameter size %d - %d",size,int(fClingParameters.size()) );
2473  for(Int_t i = 0; i < size; ++i)
2474  {
2475  TString name = TString::Format("%d",i);
2476  SetParameter(name,params[i]);
2477  }
2478  return;
2479  }
2480  fAllParametersSetted = true;
2481  std::copy(params, params+size, fClingParameters.begin() );
2482 }
2483 
2484 ////////////////////////////////////////////////////////////////////////////////
2485 /// Set a vector of parameters value.
2486 /// Order in the vector is by default the aphabetic order given to the parameters
2487 /// apart if the users has defined explicitly the parameter names
2488 
2490 {
2491  DoSetParameters(params,fNpar);
2492 }
2493 
2494 ////////////////////////////////////////////////////////////////////////////////
2495 /// Set a list of parameters.
2496 /// The order is by default the aphabetic order given to the parameters
2497 /// apart if the users has defined explicitly the parameter names
2498 
2500  Double_t p5,Double_t p6,Double_t p7,Double_t p8,
2501  Double_t p9,Double_t p10)
2502 {
2503  if(fNpar >= 1) SetParameter(0,p0);
2504  if(fNpar >= 2) SetParameter(1,p1);
2505  if(fNpar >= 3) SetParameter(2,p2);
2506  if(fNpar >= 4) SetParameter(3,p3);
2507  if(fNpar >= 5) SetParameter(4,p4);
2508  if(fNpar >= 6) SetParameter(5,p5);
2509  if(fNpar >= 7) SetParameter(6,p6);
2510  if(fNpar >= 8) SetParameter(7,p7);
2511  if(fNpar >= 9) SetParameter(8,p8);
2512  if(fNpar >= 10) SetParameter(9,p9);
2513  if(fNpar >= 11) SetParameter(10,p10);
2514 }
2515 
2516 ////////////////////////////////////////////////////////////////////////////////
2517 /// Set a parameter given a parameter index
2518 /// The parameter index is by default the aphabetic order given to the parameters
2519 /// apart if the users has defined explicitly the parameter names
2520 
2522 {
2523  if (param < 0 || param >= fNpar) return;
2524  assert(int(fClingParameters.size()) == fNpar);
2525  fClingParameters[param] = value;
2526  // TString name = TString::Format("%d",param);
2527  // SetParameter(name,value);
2528 }
2529 
2530 ////////////////////////////////////////////////////////////////////////////////
2531 void TFormula::SetParNames(const char *name0,const char *name1,const char *name2,const char *name3,
2532  const char *name4, const char *name5,const char *name6,const char *name7,
2533  const char *name8,const char *name9,const char *name10)
2534 {
2535  if(fNpar >= 1) SetParName(0,name0);
2536  if(fNpar >= 2) SetParName(1,name1);
2537  if(fNpar >= 3) SetParName(2,name2);
2538  if(fNpar >= 4) SetParName(3,name3);
2539  if(fNpar >= 5) SetParName(4,name4);
2540  if(fNpar >= 6) SetParName(5,name5);
2541  if(fNpar >= 7) SetParName(6,name6);
2542  if(fNpar >= 8) SetParName(7,name7);
2543  if(fNpar >= 9) SetParName(8,name8);
2544  if(fNpar >= 10) SetParName(9,name9);
2545  if(fNpar >= 11) SetParName(10,name10);
2546 }
2547 
2548 ////////////////////////////////////////////////////////////////////////////////
2549 void TFormula::SetParName(Int_t ipar, const char * name)
2550 {
2551 
2552  if (ipar < 0 || ipar > fNpar) {
2553  Error("SetParName","Wrong Parameter index %d ",ipar);
2554  return;
2555  }
2556  TString oldName;
2557  // find parameter with given index
2558  for ( auto &it : fParams) {
2559  if (it.second == ipar) {
2560  oldName = it.first;
2561  fParams.erase(oldName);
2562  fParams.insert(std::make_pair(name, ipar) );
2563  break;
2564  }
2565  }
2566  if (oldName.IsNull() ) {
2567  Error("SetParName","Parameter %d is not existing.",ipar);
2568  return;
2569  }
2570 
2571  //replace also parameter name in formula expression
2572  ReplaceParamName(fFormula, oldName, name);
2573 
2574 }
2575 
2576 ////////////////////////////////////////////////////////////////////////////////
2577 /// Replace in Formula expression the parameter name.
2578 
2579 void TFormula::ReplaceParamName(TString & formula, const TString & oldName, const TString & name){
2580  if (!formula.IsNull() ) {
2581  bool found = false;
2582  for(list<TFormulaFunction>::iterator it = fFuncs.begin(); it != fFuncs.end(); ++it)
2583  {
2584  if(oldName == it->GetName())
2585  {
2586  found = true;
2587  it->fName = name;
2588  break;
2589  }
2590  }
2591  if(!found)
2592  {
2593  Error("SetParName","Parameter %s is not defined.",oldName.Data());
2594  return;
2595  }
2596  // change whitespace to \\s avoid problems in parsing
2597  TString newName = name;
2598  newName.ReplaceAll(" ","\\s");
2599  TString pattern = TString::Format("[%s]",oldName.Data());
2600  TString replacement = TString::Format("[%s]",newName.Data());
2601  formula.ReplaceAll(pattern,replacement);
2602  }
2603 
2604 }
2605 
2606 ////////////////////////////////////////////////////////////////////////////////
2607 Double_t TFormula::EvalPar(const Double_t *x,const Double_t *params) const
2608 {
2609  return DoEval(x, params);
2610 }
2611 
2612 ////////////////////////////////////////////////////////////////////////////////
2613 /// Sets first 4 variables (e.g. x, y, z, t) and evaluate formula.
2614 
2616 {
2617  double xxx[4] = {x,y,z,t};
2618  return DoEval(xxx);
2619 }
2620 
2621 ////////////////////////////////////////////////////////////////////////////////
2622 /// Sets first 3 variables (e.g. x, y, z) and evaluate formula.
2623 
2625 {
2626  double xxx[3] = {x,y,z};
2627  return DoEval(xxx);
2628 }
2629 
2630 ////////////////////////////////////////////////////////////////////////////////
2631 /// Sets first 2 variables (e.g. x and y) and evaluate formula.
2632 
2634 {
2635  double xxx[2] = {x,y};
2636  return DoEval(xxx);
2637 }
2638 
2639 ////////////////////////////////////////////////////////////////////////////////
2640 /// Sets first variable (e.g. x) and evaluate formula.
2641 
2643 {
2644  //double xxx[1] = {x};
2645  double * xxx = &x;
2646  return DoEval(xxx);
2647 }
2648 
2649 ////////////////////////////////////////////////////////////////////////////////
2650 /// Evaluate formula.
2651 /// If formula is not ready to execute(missing parameters/variables),
2652 /// print these which are not known.
2653 /// If parameter has default value, and has not been setted, appropriate warning is shown.
2654 
2655 Double_t TFormula::DoEval(const double * x, const double * params) const
2656 {
2657  if(!fReadyToExecute)
2658  {
2659  Error("Eval","Formula is invalid and not ready to execute ");
2660  for(auto it = fFuncs.begin(); it != fFuncs.end(); ++it)
2661  {
2662  TFormulaFunction fun = *it;
2663  if(!fun.fFound)
2664  {
2665  printf("%s is unknown.\n",fun.GetName());
2666  }
2667  }
2668  return TMath::QuietNaN();
2669  }
2670  if (fLambdaPtr && TestBit(TFormula::kLambda)) {// case of lambda functions
2671  std::function<double(double *, double *)> & fptr = * ( (std::function<double(double *, double *)> *) fLambdaPtr);
2672  assert(x);
2673  //double * v = (x) ? const_cast<double*>(x) : const_cast<double*>(fClingVariables.data());
2674  double * v = const_cast<double*>(x);
2675  double * p = (params) ? const_cast<double*>(params) : const_cast<double*>(fClingParameters.data());
2676  return fptr(v, p);
2677  }
2678  // this is needed when reading from a file
2679  if (!fClingInitialized) {
2680  Error("Eval","Formula is invalid or not properly initialized - try calling TFormula::Compile");
2681  return TMath::QuietNaN();
2682 #ifdef EVAL_IS_NOT_CONST
2683  // need to replace in cling the name of the pointer of this object
2684  TString oldClingName = fClingName;
2685  fClingName.Replace(fClingName.Index("_0x")+1,fClingName.Length(), TString::Format("%p",this) );
2686  fClingInput.ReplaceAll(oldClingName, fClingName);
2688 #endif
2689  }
2690 
2691  Double_t result = 0;
2692  void* args[2];
2693  double * vars = (x) ? const_cast<double*>(x) : const_cast<double*>(fClingVariables.data());
2694  args[0] = &vars;
2695  if (fNpar <= 0)
2696  (*fFuncPtr)(0, 1, args, &result);
2697  else {
2698  double * pars = (params) ? const_cast<double*>(params) : const_cast<double*>(fClingParameters.data());
2699  args[1] = &pars;
2700  (*fFuncPtr)(0, 2, args, &result);
2701  }
2702  return result;
2703 }
2704 
2705 ////////////////////////////////////////////////////////////////////////////////
2706 /// Return the expression formula.
2707 ///
2708 /// - If option = "P" replace the parameter names with their values
2709 /// - If option = "CLING" return the actual expression used to build the function passed to cling
2710 /// - If option = "CLINGP" replace in the CLING expression the parameter with their values
2711 
2712 TString TFormula::GetExpFormula(Option_t *option) const
2713 {
2714  TString opt(option);
2715  if (opt.IsNull() || TestBit(TFormula::kLambda) ) return fFormula;
2716  opt.ToUpper();
2717 
2718  // if (opt.Contains("N") ) {
2719  // TString formula = fFormula;
2720  // ReplaceParName(formula, ....)
2721  // }
2722 
2723  if (opt.Contains("CLING") ) {
2724  std::string clingFunc = fClingInput.Data();
2725  std::size_t found = clingFunc.find("return");
2726  std::size_t found2 = clingFunc.rfind(";");
2727  if (found == std::string::npos || found2 == std::string::npos) {
2728  Error("GetExpFormula","Invalid Cling expression - return default formula expression");
2729  return fFormula;
2730  }
2731  TString clingFormula = fClingInput(found+7,found2-found-7);
2732  // to be implemented
2733  if (!opt.Contains("P")) return clingFormula;
2734  // replace all "p[" with "[parname"
2735  int i = 0;
2736  while (i < clingFormula.Length()-2 ) {
2737  // look for p[number
2738  if (clingFormula[i] == 'p' && clingFormula[i+1] == '[' && isdigit(clingFormula[i+2]) ) {
2739  int j = i+3;
2740  while ( isdigit(clingFormula[j]) ) { j++;}
2741  if (clingFormula[j] != ']') {
2742  Error("GetExpFormula","Parameters not found - invalid expression - return default cling formula");
2743  return clingFormula;
2744  }
2745  TString parNumbName = clingFormula(i+2,j-i-2);
2746  int parNumber = parNumbName.Atoi();
2747  assert(parNumber < fNpar);
2748  TString replacement = TString::Format("%f",GetParameter(parNumber));
2749  clingFormula.Replace(i,j-i+1, replacement );
2750  i += replacement.Length();
2751  }
2752  i++;
2753  }
2754  return clingFormula;
2755  }
2756  if (opt.Contains("P") ) {
2757  // replace parameter names with their values
2758  TString expFormula = fFormula;
2759  int i = 0;
2760  while (i < expFormula.Length()-2 ) {
2761  // look for [parName]
2762  if (expFormula[i] == '[') {
2763  int j = i+1;
2764  while ( expFormula[j] != ']' ) { j++;}
2765  if (expFormula[j] != ']') {
2766  Error("GetExpFormula","Parameter names not found - invalid expression - return default formula");
2767  return expFormula;
2768  }
2769  TString parName = expFormula(i+1,j-i-1);
2770  TString replacement = TString::Format("%g",GetParameter(parName));
2771  expFormula.Replace(i,j-i+1, replacement );
2772  i += replacement.Length();
2773  }
2774  i++;
2775  }
2776  return expFormula;
2777  }
2778  Warning("GetExpFormula","Invalid option - return defult formula expression");
2779  return fFormula;
2780 }
2781 
2782 ////////////////////////////////////////////////////////////////////////////////
2783 /// Print the formula and its attributes.
2784 
2785 void TFormula::Print(Option_t *option) const
2786 {
2787  printf(" %20s : %s Ndim= %d, Npar= %d, Number= %d \n",GetName(),GetTitle(), fNdim,fNpar,fNumber);
2788  printf(" Formula expression: \n");
2789  printf("\t%s \n",fFormula.Data() );
2790  TString opt(option);
2791  opt.ToUpper();
2792  // do an evaluation as a cross-check
2793  //if (fReadyToExecute) Eval();
2794 
2795  if (opt.Contains("V") ) {
2796  if (fNdim > 0 && !TestBit(TFormula::kLambda)) {
2797  printf("List of Variables: \n");
2798  assert(int(fClingVariables.size()) >= fNdim);
2799  for ( int ivar = 0; ivar < fNdim ; ++ivar) {
2800  printf("Var%4d %20s = %10f \n",ivar,GetVarName(ivar).Data(), fClingVariables[ivar]);
2801  }
2802  }
2803  if (fNpar > 0) {
2804  printf("List of Parameters: \n");
2805  if ( int(fClingParameters.size()) < fNpar)
2806  Error("Print","Number of stored parameters in vector %zu in map %zu is different than fNpar %d",fClingParameters.size(), fParams.size(), fNpar);
2807  assert(int(fClingParameters.size()) >= fNpar);
2808  // print with order passed to Cling function
2809  for ( int ipar = 0; ipar < fNpar ; ++ipar) {
2810  printf("Par%4d %20s = %10f \n",ipar,GetParName(ipar), fClingParameters[ipar] );
2811  }
2812  }
2813  printf("Expression passed to Cling:\n");
2814  printf("\t%s\n",fClingInput.Data() );
2815  }
2816  if(!fReadyToExecute)
2817  {
2818  Warning("Print","Formula is not ready to execute. Missing parameters/variables");
2819  for(list<TFormulaFunction>::const_iterator it = fFuncs.begin(); it != fFuncs.end(); ++it)
2820  {
2821  TFormulaFunction fun = *it;
2822  if(!fun.fFound)
2823  {
2824  printf("%s is unknown.\n",fun.GetName());
2825  }
2826  }
2827  }
2829  {
2830  // we can skip this
2831  // Info("Print","Not all parameters are setted.");
2832  // for(map<TString,TFormulaVariable>::const_iterator it = fParams.begin(); it != fParams.end(); ++it)
2833  // {
2834  // pair<TString,TFormulaVariable> param = *it;
2835  // if(!param.second.fFound)
2836  // {
2837  // printf("%s has default value %lf\n",param.first.Data(),param.second.GetInitialValue());
2838  // }
2839  // }
2840 
2841  }
2842 
2843 
2844 }
2845 
2846 ////////////////////////////////////////////////////////////////////////////////
2847 /// Stream a class object.
2848 
2849 void TFormula::Streamer(TBuffer &b)
2850 {
2851  if (b.IsReading() ) {
2852  UInt_t R__s, R__c;
2853  Version_t v = b.ReadVersion(&R__s, &R__c);
2854  //std::cout << "version " << v << std::endl;
2855  if (v <= 8 && v > 3 && v != 6) {
2856  // old TFormula class
2857  ROOT::v5::TFormula * fold = new ROOT::v5::TFormula();
2858  // read old TFormula class
2859  fold->Streamer(b, v, R__s, R__c, TFormula::Class());
2860  //std::cout << "read old tformula class " << std::endl;
2861  TFormula fnew(fold->GetName(), fold->GetExpFormula() );
2862 
2863  *this = fnew;
2864 
2865 // printf("copying content in a new TFormula \n");
2866  SetParameters(fold->GetParameters() );
2867  if (!fReadyToExecute ) {
2868  Error("Streamer","Old formula read from file is NOT valid");
2869  Print("v");
2870  }
2871  delete fold;
2872  return;
2873  }
2874  else if (v > 8) {
2875  // new TFormula class
2876  b.ReadClassBuffer(TFormula::Class(), this, v, R__s, R__c);
2877 
2878  //std::cout << "reading npar = " << GetNpar() << std::endl;
2879 
2880  // initialize the formula
2881  // need to set size of fClingVariables which is transient
2882  //fClingVariables.resize(fNdim);
2883 
2884  // case of formula contains only parameters
2885  if (fFormula.IsNull() ) return;
2886 
2887 
2888  // store parameter values, names and order
2889  std::vector<double> parValues = fClingParameters;
2890  auto paramMap = fParams;
2891  fNpar = fParams.size();
2892 
2893  if (!TestBit(TFormula::kLambda) ) {
2894 
2895  //std::cout << "Streamer::Reading preprocess the formula " << fFormula << " ndim = " << fNdim << " npar = " << fNpar << std::endl;
2896  // for ( auto &p : fParams)
2897  // std::cout << "parameter " << p.first << " index " << p.second << std::endl;
2898 
2899  fClingParameters.clear(); // need to be reset before re-initializing it
2900 
2901  FillDefaults();
2902 
2903 
2905 
2906  //std::cout << "Streamer::after pre-process the formula " << fFormula << " ndim = " << fNdim << " npar = " << fNpar << std::endl;
2907 
2909 
2910  //std::cout << "Streamer::after prepared " << fClingInput << " ndim = " << fNdim << " npar = " << fNpar << std::endl;
2911 
2912 
2913  // restore parameter values
2914  if (fNpar != (int) parValues.size() ) {
2915  Error("Streamer","number of parameters computed (%d) is not same as the stored parameters (%d)",fNpar,int(parValues.size()) );
2916  Print("v");
2917  }
2918  }
2919  else {
2920  // case of lamda expressions
2921  bool ret = InitLambdaExpression(fFormula);
2922  if (ret) {
2923  fReadyToExecute = true;
2924  fClingInitialized = true;
2925  }
2926  }
2927  assert(fNpar == (int) parValues.size() );
2928  std::copy( parValues.begin(), parValues.end(), fClingParameters.begin() );
2929  // restore parameter names and order
2930  if (fParams.size() != paramMap.size() ) {
2931  Warning("Streamer","number of parameters list found (%zu) is not same as the stored one (%zu) - use re-created list",fParams.size(),paramMap.size()) ;
2932  //Print("v");
2933  }
2934  else
2935  //assert(fParams.size() == paramMap.size() );
2936  fParams = paramMap;
2937 
2938  // input formula into Cling
2939  // need to replace in cling the name of the pointer of this object
2940  // TString oldClingName = fClingName;
2941  // fClingName.Replace(fClingName.Index("_0x")+1,fClingName.Length(), TString::Format("%p",this) );
2942  // fClingInput.ReplaceAll(oldClingName, fClingName);
2943  // InputFormulaIntoCling();
2944 
2945 
2946 
2947  if (!TestBit(kNotGlobal)) {
2949  gROOT->GetListOfFunctions()->Add(this);
2950  }
2951  if (!fReadyToExecute ) {
2952  Error("Streamer","Formula read from file is NOT ready to execute");
2953  Print("v");
2954  }
2955  //std::cout << "reading 2 npar = " << GetNpar() << std::endl;
2956 
2957  return;
2958  }
2959  else {
2960  Error("Streamer","Reading version %d is not supported",v);
2961  return;
2962  }
2963  }
2964  else {
2965  // case of writing
2967  //std::cout << "writing npar = " << GetNpar() << std::endl;
2968  }
2969 }
constexpr Double_t H()
Definition: TMath.h:140
virtual void Clear(Option_t *option="")
Clear the formula setting expression to empty and reset the variables and parameters containers...
Definition: TFormula.cxx:621
static Bool_t IsBracket(const char c)
Definition: TFormula.cxx:153
virtual const char * GetName() const
Returns name of object.
Definition: TNamed.h:47
Bool_t IsReading() const
Definition: TBuffer.h:81
bool operator()(const TString &a, const TString &b) const
Definition: TFormula.cxx:233
void Streamer(TBuffer &b, const TClass *onfile_class)
Stream a class object.
void HandleExponentiation(TString &formula)
Handling exponentiation Can handle multiple carets, eg.
Definition: TFormula.cxx:1199
std::map< TString, Double_t > fConsts
Definition: TFormula.h:115
constexpr Double_t K()
Definition: TMath.h:178
virtual Int_t WriteClassBuffer(const TClass *cl, void *pointer)=0
virtual void Info(const char *method, const char *msgfmt,...) const
Issue info message.
Definition: TObject.cxx:847
Double_t Eval(Double_t x) const
Sets first variable (e.g. x) and evaluate formula.
Definition: TFormula.cxx:2642
Int_t fNpar
Definition: TFormula.h:119
Double_t Floor(Double_t x)
Definition: TMath.h:600
virtual TFormula * GetFormula()
Definition: TF1.h:407
std::vector< TObject * > fLinearParts
Definition: TFormula.h:121
static double p3(double t, double a, double b, double c, double d)
short Version_t
Definition: RtypesCore.h:61
virtual CallFuncIFacePtr_t CallFunc_IFacePtr(CallFunc_t *) const
Definition: TInterpreter.h:272
static const TString gNamePrefix
Definition: TFormula.cxx:134
TString GetVarName(Int_t ivar) const
Returns variable name given its position in the array.
Definition: TFormula.cxx:2218
void AddVariables(const TString *vars, const Int_t size)
Adds multiple variables.
Definition: TFormula.cxx:2098
TMethodCall * fMethod
Definition: TFormula.h:95
constexpr Double_t Sqrt2()
Definition: TMath.h:68
const char Option_t
Definition: RtypesCore.h:62
void DoSetParameters(const Double_t *p, Int_t size)
Definition: TFormula.cxx:2467
void HandleLinear(TString &formula)
Handle linear functions defined with the operator ++.
Definition: TFormula.cxx:1296
void HandlePolN(TString &formula)
Handling polN If before &#39;pol&#39; exist any name, this name will be treated as variable used in polynomia...
Definition: TFormula.cxx:828
Double_t QuietNaN()
Definition: TMath.h:784
const Ssiz_t kNPOS
Definition: RtypesCore.h:115
TString fName
Definition: TFormula.h:29
CallFunc_t * GetCallFunc() const
Definition: TMethodCall.h:93
Bool_t TestBit(UInt_t f) const
Definition: TObject.h:159
void SetParameters(const Double_t *params)
Set a vector of parameters value.
Definition: TFormula.cxx:2489
virtual void SetName(const char *name)
Set the name of the TNamed.
Definition: TNamed.cxx:131
virtual void Copy(TObject &f1) const
Copy this to obj.
Definition: TFormula.cxx:544
void SetParNames(const char *name0="p0", const char *name1="p1", const char *name2="p2", const char *name3="p3", const char *name4="p4", const char *name5="p5", const char *name6="p6", const char *name7="p7", const char *name8="p8", const char *name9="p9", const char *name10="p10")
Definition: TFormula.cxx:2531
void SetVariables(const std::pair< TString, Double_t > *vars, const Int_t size)
Sets multiple variables.
Definition: TFormula.cxx:2170
Buffer base class used for serializing objects.
Definition: TBuffer.h:40
Helper class for TFormula.
Definition: TFormula.h:26
#define gROOT
Definition: TROOT.h:375
Int_t GetNpar() const
Definition: TFormula.h:170
void AddVariable(const TString &name, Double_t value=0)
Adds variable to known variables, and reprocess formula.
Definition: TFormula.cxx:2059
void InputFormulaIntoCling()
pointer to the lambda function
Definition: TFormula.cxx:713
void DoAddParameter(const TString &name, Double_t value, bool processFormula)
Adds parameter to known parameters.
Definition: TFormula.cxx:2250
int Int_t
Definition: RtypesCore.h:41
bool Bool_t
Definition: RtypesCore.h:59
TArc * a
Definition: textangle.C:12
R__EXTERN TVirtualMutex * gROOTMutex
Definition: TROOT.h:57
const TList * GetListOfAllPublicMethods(Bool_t load=kTRUE)
Returns a list of all public methods of this class and its base classes.
Definition: TClass.cxx:3644
Bool_t fFound
Definition: TFormula.h:32
Bool_t fAllParametersSetted
transient to force re-initialization
Definition: TFormula.h:94
#define gInterpreter
Definition: TInterpreter.h:499
STL namespace.
void SetPredefinedParamNames()
Set parameter names only in case of pre-defined functions.
Definition: TFormula.cxx:1971
void SetBit(UInt_t f, Bool_t set)
Set or unset the user status bits as specified in f.
Definition: TObject.cxx:687
Bool_t PrepareFormula(TString &formula)
prepare the formula to be executed normally is called with fFormula
Definition: TFormula.cxx:1371
Int_t GetNumber() const
Definition: TFormula.h:171
constexpr Double_t Ln10()
Definition: TMath.h:80
void ReplaceParamName(TString &formula, const TString &oldname, const TString &name)
Replace in Formula expression the parameter name.
Definition: TFormula.cxx:2579
virtual ~TFormula()
Definition: TFormula.cxx:292
Double_t x[n]
Definition: legend1.C:17
static TString Format(const char *fmt,...)
Static method which formats a string using a printf style format descriptor and return a TString...
Definition: TString.cxx:2345
Bool_t IsFuncCall() const
Definition: TFormula.h:37
void Class()
Definition: Class.C:29
Int_t fArrayPos
Definition: TFormula.h:64
The TNamed class is the base class for all named ROOT classes.
Definition: TNamed.h:29
void SetParameter(const char *name, Double_t value)
Sets parameter value.
Definition: TFormula.cxx:2405
Double_t Log10(Double_t x)
Definition: TMath.h:652
virtual Bool_t Declare(const char *code)=0
static double p2(double t, double a, double b, double c)
void(* Generic_t)(void *, int, void **, void *)
Definition: TInterpreter.h:74
TInterpreter::CallFuncIFacePtr_t::Generic_t fFuncPtr
unique name passed to Cling to define the function ( double clingName(double*x, double*p) ) ...
Definition: TFormula.h:98
static Bool_t IsScientificNotation(const TString &formula, int ipos)
Definition: TFormula.cxx:176
TFormula & operator=(const TFormula &rhs)
= operator.
Definition: TFormula.cxx:459
static bool IsReservedName(const char *name)
Definition: TFormula.cxx:283
const TObject * GetLinearPart(Int_t i) const
Return linear part.
Definition: TFormula.cxx:2043
std::vector< std::vector< double > > Data
static Bool_t IsAParameterName(const TString &formula, int ipos)
Definition: TFormula.cxx:211
Double_t fValue
Definition: TFormula.h:63
void * fLambdaPtr
function pointer
Definition: TFormula.h:99
Method or function calling interface.
Definition: TMethodCall.h:37
void SetParName(Int_t ipar, const char *name)
Definition: TFormula.cxx:2549
constexpr Double_t Pi()
Definition: TMath.h:40
static std::unordered_map< std::string, void * > gClingFunctions
Definition: TFormula.cxx:138
Double_t Infinity()
Definition: TMath.h:797
Int_t GetNdim() const
Definition: TFormula.h:169
A doubly linked list.
Definition: TList.h:43
Int_t fNdim
Definition: TFormula.h:118
Double_t GetVariable(const char *name) const
Returns variable value.
Definition: TFormula.cxx:2190
Double_t EvalPar(const Double_t *x, const Double_t *params=0) const
Definition: TFormula.cxx:2607
virtual void Delete(Option_t *option="")
Delete this object.
Definition: TObject.cxx:176
static Bool_t IsHexadecimal(const TString &formula, int ipos)
Definition: TFormula.cxx:188
Int_t GetNargs() const
Number of function arguments.
Definition: TFunction.cxx:164
constexpr Double_t G()
Definition: TMath.h:106
virtual TObject * FindObject(const char *name) const
Must be redefined in derived classes.
Definition: TObject.cxx:328
Int_t GetNargs() const
Definition: TFormula.h:36
std::vector< Double_t > fClingParameters
cached variables
Definition: TFormula.h:91
void FillDefaults()
Fill structures with default variables, constants and function shortcuts.
Definition: TFormula.cxx:726
void HandleParametrizedFunctions(TString &formula)
Handling parametrized functions Function can be normalized, and have different variable then x...
Definition: TFormula.cxx:940
SVector< double, 2 > v
Definition: Dict.h:5
static Bool_t IsOperator(const char c)
Definition: TFormula.cxx:141
void SetName(const char *name)
Set the name of the formula.
Definition: TFormula.cxx:2142
The Formula class.
Definition: TFormula.h:83
const char * GetName() const
Definition: TFormula.h:34
unsigned int UInt_t
Definition: RtypesCore.h:42
TMarker * m
Definition: textangle.C:8
static Bool_t IsDefaultVariableName(const TString &name)
Definition: TFormula.cxx:170
virtual void Error(const char *method, const char *msgfmt,...) const
Issue error message.
Definition: TObject.cxx:873
static Bool_t IsFunctionNameChar(const char c)
Definition: TFormula.cxx:164
The ROOT global object gROOT contains a list of all defined classes.
Definition: TClass.h:71
const std::string sname
Definition: testIO.cxx:45
static double p1(double t, double a, double b)
Int_t GetParNumber(const char *name) const
Return parameter index given a name (return -1 for not existing parameters) non need to print an erro...
Definition: TFormula.cxx:2331
Bool_t PrepareEvalMethod()
Sets TMethodCall to function inside Cling environment.
Definition: TFormula.cxx:658
TString fFormula
Definition: TFormula.h:117
Bool_t IsValid() const
Return true if the method call has been properly initialized and is usable.
constexpr Double_t E()
Definition: TMath.h:74
void InitWithPrototype(TClass *cl, const char *method, const char *proto, Bool_t objectIsConst=kFALSE, ROOT::EFunctionMatchMode mode=ROOT::kConversionMatch)
Initialize the method invocation environment.
#define R__LOCKGUARD2(mutex)
TString fClingName
pointer to methocall
Definition: TFormula.h:96
std::vector< Double_t > fClingVariables
input function passed to Cling
Definition: TFormula.h:90
virtual Int_t ReadClassBuffer(const TClass *cl, void *pointer, const TClass *onfile_class=0)=0
int Ssiz_t
Definition: RtypesCore.h:63
std::map< TString, TString > fFunctionsShortcuts
Definition: TFormula.h:116
#define ClassImp(name)
Definition: Rtypes.h:336
double f(double x)
double Double_t
Definition: RtypesCore.h:55
The FORMULA class (ROOT version 5)
Definition: TFormula.h:65
Double_t y[n]
Definition: legend1.C:17
const char * GetParName(Int_t ipar) const
Return parameter name given by integer.
Definition: TFormula.cxx:2370
Int_t fNumber
Definition: TFormula.h:120
Int_t Compile(const char *expression="")
Compile the given expression with Cling backward compatibility method to be used in combination with ...
Definition: TFormula.cxx:510
constexpr Double_t C()
Definition: TMath.h:92
static TClass * GetClass(const char *name, Bool_t load=kTRUE, Bool_t silent=kFALSE)
Static method returning pointer to TClass of the specified class name.
Definition: TClass.cxx:2885
Double_t DoEval(const Double_t *x, const Double_t *p=nullptr) const
Evaluate formula.
Definition: TFormula.cxx:2655
void SetVariable(const TString &name, Double_t value)
Sets variable value.
Definition: TFormula.cxx:2234
Mother of all ROOT objects.
Definition: TObject.h:37
constexpr Double_t EulerGamma()
Definition: TMath.h:238
Bool_t InitLambdaExpression(const char *formula)
Definition: TFormula.cxx:469
TString fName
Definition: TFormula.h:62
you should not use this method at all Int_t Int_t z
Definition: TRolke.cxx:630
Global functions class (global functions are obtained from CINT).
Definition: TFunction.h:28
virtual void Copy(TObject &named) const
Copy this to obj.
Definition: TNamed.cxx:85
TString GetExpFormula(Option_t *option="") const
Return the expression formula.
Definition: TFormula.cxx:2712
Short_t Max(Short_t a, Short_t b)
Definition: TMathBase.h:200
1-Dim function class
Definition: TF1.h:150
void ProcessFormula(TString &formula)
Iterates through funtors in fFuncs and performs the appropriate action.
Definition: TFormula.cxx:1637
Each ROOT class (see TClass) has a linked list of methods.
Definition: TMethod.h:38
TF1 * f1
Definition: legend1.C:11
you should not use this method at all Int_t Int_t Double_t Double_t Double_t Int_t Double_t Double_t Double_t Double_t b
Definition: TRolke.cxx:630
virtual Double_t * GetParameters() const
Definition: TFormula.h:243
constexpr Double_t LogE()
Definition: TMath.h:86
R__EXTERN Int_t gDebug
Definition: Rtypes.h:83
Double_t GetParameter(const char *name) const
Returns parameter value given by string.
Definition: TFormula.cxx:2344
constexpr Double_t Sigma()
Definition: TMath.h:192
double result[121]
Bool_t fReadyToExecute
Definition: TFormula.h:92
Definition: first.py:1
virtual TString GetExpFormula(Option_t *option="") const
Reconstruct the formula expression from the internal TFormula member variables.
R__EXTERN TInterpreter * gCling
Definition: TInterpreter.h:501
Int_t GetVarNumber(const char *name) const
Returns variable number (positon in array) given its name.
Definition: TFormula.cxx:2204
void ExtractFunctors(TString &formula)
Extracts functors from formula, and put them in fFuncs.
Definition: TFormula.cxx:1410
TString fClingInput
Definition: TFormula.h:89
void Print(Option_t *option="") const
Print the formula and its attributes.
Definition: TFormula.cxx:2785
Bool_t IsValid() const
Definition: TFormula.h:181
constexpr Double_t R()
Definition: TMath.h:213
std::map< TString, TFormulaVariable > fVars
Definition: TFormula.h:113
const Int_t n
Definition: legend1.C:16
Int_t GetNargsOpt() const
Number of function optional (default) arguments.
Definition: TFunction.cxx:174
Bool_t fClingInitialized
trasient to force initialization
Definition: TFormula.h:93
Double_t * GetParameters() const
Definition: TFormula.cxx:2384
std::map< TString, Int_t, TFormulaParamOrder > fParams
list of variable names
Definition: TFormula.h:114
void variables(TString dataset, TString fin="TMVA.root", TString dirName="InputVariables_Id", TString title="TMVA Input Variables", Bool_t isRegression=kFALSE, Bool_t useTMVAStyle=kTRUE)
std::list< TFormulaFunction > fFuncs
Definition: TFormula.h:112
Another helper class for TFormula.
Definition: TFormula.h:59
virtual void Warning(const char *method, const char *msgfmt,...) const
Issue warning message.
Definition: TObject.cxx:859
virtual Version_t ReadVersion(UInt_t *start=0, UInt_t *bcnt=0, const TClass *cl=0)=0
virtual const char * GetTitle() const
Returns title of object.
Definition: TNamed.h:48
void PreProcessFormula(TString &formula)
Preprocessing of formula Replace all ** by ^, and removes spaces.
Definition: TFormula.cxx:1351