Logo ROOT  

Reference Guide
 
Loading...
Searching...
No Matches
TFormula_v5.cxx
Go to the documentation of this file.
1// @(#)root/hist:$Id$
2// Author: Nicolas Brun 19/08/95
3
4/*************************************************************************
5 * Copyright (C) 1995-2000, Rene Brun and Fons Rademakers. *
6 * All rights reserved. *
7 * *
8 * For the licensing terms see $ROOTSYS/LICENSE. *
9 * For the list of contributors see $ROOTSYS/README/CREDITS. *
10 *************************************************************************/
11
12#include <cmath>
13
14#include "TROOT.h"
15#include "TClass.h"
16#include "TBuffer.h"
17#include "v5/TFormula.h"
18#include "TMath.h"
19#include "TRandom.h"
20#include "TFunction.h"
21#include "TMethodCall.h"
22#include "TObjString.h"
23#include "TError.h"
24#include "v5/TFormulaPrimitive.h"
25#include "TInterpreter.h"
26#include "TVirtualMutex.h"
27#include "strlcpy.h"
28#include "snprintf.h"
29
30#ifdef WIN32
31#pragma optimize("",off)
32#endif
33
34static Int_t gMAXOP=1000,gMAXPAR=1000,gMAXCONST=1000;
35const Int_t gMAXSTRINGFOUND = 10;
36const UInt_t kOptimizationError = BIT(19);
37
38ClassImp(ROOT::v5::TFormula);
39
40namespace ROOT {
41
42 namespace v5 {
43
44/** \class TFormula TFormula.h "inc/v5/TFormula.h"
45 \ingroup Hist
46The FORMULA class (ROOT version 5)
47
48 Example of valid expressions:
49
50 - `sin(x)/x`
51 - `[0]*sin(x) + [1]*exp(-[2]*x)`
52 - `x + y**2`
53 - `x^2 + y^2`
54 - `[0]*pow([1],4)`
55 - `2*pi*sqrt(x/y)`
56 - `gaus(0)*expo(3) + ypol3(5)*x`
57 - `gausn(0)*expo(3) + ypol3(5)*x`
58
59 In the last example above:
60
61 gaus(0) is a substitute for `[0]*exp(-0.5*((x-[1])/[2])**2)`
62 and (0) means start numbering parameters at 0
63
64 gausn(0) is a substitute for `[0]*exp(-0.5*((x-[1])/[2])**2)/(sqrt(2*pi)*[2]))`
65 and (0) means start numbering parameters at 0
66
67 expo(3) is a substitute for `exp([3]+[4]*x)`
68
69 pol3(5) is a substitute for `par[5]+par[6]*x+par[7]*x**2+par[8]*x**3`
70 (here Pol3 stands for Polynomial of degree 3)
71
72 TMath functions can be part of the expression, eg:
73
74 - `TMath::Landau(x)*sin(x)`
75 - `TMath::Erf(x)`
76
77 Comparisons operators are also supported (&&, ||, ==, <=, >=, !)
78 Examples:
79
80 sin(x*(x<0.5 || x>1))
81
82 If the result of a comparison is TRUE, the result is 1, otherwise 0.
83
84 Already predefined names can be given. For example, if the formula
85
86 TFormula old(sin(x*(x<0.5 || x>1))) one can assign a name to the formula. By default
87 the name of the object = title = formula itself.
88 old.SetName("old").
89 then, old can be reused in a new expression.
90 TFormula new("x*old") is equivalent to:
91 TFormula new("x*sin(x*(x<0.5 || x>1))")
92
93 Up to 4 dimensions are supported (indicated by x, y, z, t)
94 An expression may have 0 parameters or a list of parameters
95 indicated by the sequence [par_number]
96
97 A graph showing the logic to compile and analyze a formula
98 is shown in TFormula::Compile and TFormula::Analyze.
99 Once a formula has been compiled, it can be evaluated for a given
100 set of parameters. see graph in TFormula::EvalPar.
101
102 This class is the base class for the function classes TF1,TF2 and TF3.
103 It is also used by the ntuple selection mechanism TNtupleFormula.
104
105 In version 7 of TFormula, the usage of fOper has been changed
106 to improve the performance of TFormula::EvalPar.
107 Conceptually, fOper was changed from a simple array of Int_t
108 to an array of composite values.
109 For example a 'ylandau(5)' operation used to be encoded as 4105;
110 it is now encoded as (klandau >> kTFOperShift) + 5
111 Any class inheriting from TFormula and using directly fOper (which
112 is now a private data member), needs to be updated to take this
113 in consideration. The member functions recommended to set and
114 access fOper are: SetAction, GetAction, GetActionParam
115 For more performant access to the information, see the implementation
116 TFormula::EvalPar
117
118 ### CHANGING DEFAULT SETTINGS
119
120 When creating complex formula , it may be necessary to increase
121 some default parameters. see static function TFormula::SetMaxima
122
123 ### WHY TFormula CANNOT ACCEPT A CLASS MEMBER FUNCTION ?
124
125 This is a frequently asked question.
126 C++ is a strongly typed language. There is no way for TFormula (without
127 recompiling this class) to know about all possible user defined data types.
128 This also apply to the case of a static class function.
129 Because TMath is a special and frequent case, TFormula is aware
130 of all TMath functions.
131*/
132
133////////////////////////////////////////////////////////////////////////////////
134/// Formula default constructor.
135
136TFormula::TFormula(): TNamed()
137{
138 fNdim = 0;
139 fNpar = 0;
140 fNoper = 0;
141 fNconst = 0;
142 fNumber = 0;
143 fExpr = 0;
144 fOper = 0;
145 fConst = 0;
146 fParams = 0;
147 fNstring= 0;
148 fNames = 0;
149 fNval = 0;
150 //
151 //MI change
152 fNOperOptimized = 0;
153 fExprOptimized = 0;
154 fOperOptimized = 0;
155 fOperOffset = 0;
156 fPredefined = 0;
157 fOptimal = (ROOT::v5::TFormulaPrimitive::TFuncG)&TFormula::EvalParOld;
158}
159
160////////////////////////////////////////////////////////////////////////////////
161/// Normal Formula constructor.
162
163TFormula::TFormula(const char *name,const char *expression) :
164 TNamed(name,expression)
165{
166 fNdim = 0;
167 fNpar = 0;
168 fNoper = 0;
169 fNconst = 0;
170 fNumber = 0;
171 fExpr = 0;
172 fOper = 0;
173 fConst = 0;
174 fParams = 0;
175 fNstring= 0;
176 fNames = 0;
177 fNval = 0;
178 //
179 //MI change
180 fNOperOptimized = 0;
181 fExprOptimized = 0;
182 fOperOptimized = 0;
183 fOperOffset = 0;
184 fPredefined = 0;
185 fOptimal = (ROOT::v5::TFormulaPrimitive::TFuncG)&TFormula::EvalParOld;
186
187 if (!expression || !*expression) {
188 Error("TFormula", "expression may not be 0 or have 0 length");
189 return;
190 }
191
192 //eliminate blanks in expression
193 Int_t i,j,nch;
194 nch = strlen(expression);
195 char *expr = new char[nch+1];
196 j = 0;
197 for (i=0;i<nch;i++) {
198 if (expression[i] == ' ') continue;
199 if (i > 0 && (expression[i] == '*') && (expression[i-1] == '*')) {
200 expr[j-1] = '^';
201 continue;
202 }
203 expr[j] = expression[i]; j++;
204 }
205 expr[j] = 0;
206 Bool_t gausNorm = kFALSE;
207 Bool_t landauNorm = kFALSE;
208 Bool_t linear = kFALSE;
209
210 if (j) {
211 TString chaine = expr;
212 //special case for functions for linear fitting
213 if (chaine.Contains("++"))
214 linear = kTRUE;
215 // special case for normalized gaus
216 if (chaine.Contains("gausn")) {
217 gausNorm = kTRUE;
218 TString tmp = chaine;
219 tmp.ReplaceAll("gausn","");
220 tmp.ReplaceAll("landaun","");
221 if ( tmp.Contains("gaus") )
222 Warning("TFormula","Cannot use both gaus and gausn - gaus will be treated as gausn");
223 if ( tmp.Contains("landau") )
224 Warning("TFormula","Cannot use both gausn and landau - landau will be treated as landaun");
225 }
226 // special case for normalized landau
227 if (chaine.Contains("landaun")) {
228 landauNorm = kTRUE;
229 TString tmp = chaine;
230 tmp.ReplaceAll("landaun","");
231 tmp.ReplaceAll("gausn","");
232 if ( tmp.Contains("gaus") ) {
233 Warning("TFormula","Cannot use both gaus and landaun - gaus will be treated as gausn");
234 }
235 if ( tmp.Contains("landau") )
236 Warning("TFormula","Cannot use both landau and landaun - landau will be treated as landaun");
237 }
238 // need to the replacement here for the error message before
239 if (gausNorm)
240 chaine.ReplaceAll("gausn","gaus");
241 if (landauNorm)
242 chaine.ReplaceAll("landaun","landau");
243
244 SetTitle(chaine.Data());
245 }
246 delete [] expr;
247
248 if (linear) SetBit(kLinear);
249
250 if (Compile()) return;
251
252 if (gausNorm) SetBit(kNormalized);
253 if (landauNorm) SetBit(kNormalized);
254
255 // Store formula in linked list of formula in ROOT
256
257
258 if (strcmp(name,"x")==0 || strcmp(name,"y")==0 ||
259 strcmp(name,"z")==0 || strcmp(name,"t")==0 )
260 {
261 Error("TFormula","The name \'%s\' is reserved as a TFormula variable name.\n"
262 "\tThis function will not be registered in the list of functions",name);
263 } else {
264 R__LOCKGUARD(gROOTMutex);
265 TFormula *old = (TFormula*)gROOT->GetListOfFunctions()->FindObject(name);
266 if (old) {
267 gROOT->GetListOfFunctions()->Remove(old);
268 }
269 gROOT->GetListOfFunctions()->Add(this);
270 }
271}
272
273////////////////////////////////////////////////////////////////////////////////
274/// Default constructor.
275
276TFormula::TFormula(const TFormula &formula) : TNamed()
277{
278 fNdim = 0;
279 fNpar = 0;
280 fNoper = 0;
281 fNconst = 0;
282 fNumber = 0;
283 fExpr = 0;
284 fOper = 0;
285 fConst = 0;
286 fParams = 0;
287 fNstring= 0;
288 fNames = 0;
289 fNval = 0;
290 fNOperOptimized = 0;
291 fPredefined = 0;
292 fOperOffset = 0;
293 fExprOptimized = 0;
294 fOperOptimized = 0;
295 fOptimal = (ROOT::v5::TFormulaPrimitive::TFuncG)&TFormula::EvalParOld;
296
297 ((TFormula&)formula).TFormula::Copy(*this);
298}
299
300////////////////////////////////////////////////////////////////////////////////
301/// Operator =
302
303TFormula& TFormula::operator=(const TFormula &rhs)
304{
305 if (this != &rhs) {
306 rhs.Copy(*this);
307 }
308 return *this;
309}
310
311////////////////////////////////////////////////////////////////////////////////
312/// Formula default destructor.
313
314TFormula::~TFormula()
315{
316 if (gROOT) {
317 R__LOCKGUARD(gROOTMutex);
318 gROOT->GetListOfFunctions()->Remove(this);
319 }
320
321 ClearFormula();
322}
323
324////////////////////////////////////////////////////////////////////////////////
325/// Check if the chain as function call.
326///
327/// If you overload this member function, you also HAVE TO
328/// never call the constructor:
329///
330/// ~~~ {.cpp}
331/// TFormula::TFormula(const char *name,const char *expression)
332/// ~~~
333///
334/// and write your own constructor
335///
336/// ~~~ {.cpp}
337/// MyClass::MyClass(const char *name,const char *expression) : TFormula()
338/// ~~~
339///
340/// which has to call the TFormula default constructor and whose implementation
341/// should be similar to the implementation of the normal TFormula constructor
342///
343/// This is necessary because the normal TFormula constructor call indirectly
344/// the virtual member functions Analyze, DefaultString, DefaultValue
345/// and DefaultVariable.
346
347Bool_t TFormula::AnalyzeFunction(TString &chaine, Int_t &err, Int_t offset)
348{
349 int i;
350
351 // We have to decompose the chain is 3 potential components:
352 // namespace::functionName( args )
353
354 Ssiz_t argStart = chaine.First('(');
355 if (argStart<0) return false;
356
357 TString functionName = chaine(0,argStart);
358
359 // This does not support template yet (where the scope operator might be in
360 // one of the template arguments
361 Ssiz_t scopeEnd = functionName.Last(':');
362 TString spaceName;
363 if (scopeEnd>0 && functionName[scopeEnd-1]==':') {
364 spaceName = functionName(0,scopeEnd-1);
365 functionName.Remove(0,scopeEnd+1);
366 }
367
368 // Now we need to count and decompose the actual arguments, we could also check the type
369 // of the arguments
370 if (chaine[chaine.Length()-1] != ')') {
371 Error("AnalyzeFunction","We thought we had a function but we dont (in %s)\n",chaine.Data());
372 }
373
374 TString args = chaine(argStart+1,chaine.Length()-2-argStart);
375 TObjArray argArr;
376 argArr.SetOwner(kTRUE);
377 //fprintf(stderr,"args are '%s'\n",args.Data());
378
379 Bool_t inString = false;
380 int paran = 0;
381 int brack = 0;
382 int prevComma = 0;
383 int nargs = 0;
384 for(i=0; i<args.Length(); i++) {
385 if (args[i]=='"') inString = !inString;
386 if (inString) continue;
387
388 Bool_t foundArg = false;
389 switch(args[i]) {
390
391 case '(': paran++; break;
392 case ')': paran--; break;
393 case '[': brack++; break;
394 case ']': brack--; break;
395
396 case ',': if (paran==0 && brack==0) { foundArg = true; } break;
397 }
398 if ((i+1)==args.Length()) {
399 foundArg = true; i++;
400 }
401 if (foundArg) {
402 TString arg = args(prevComma,i-prevComma);
403
404 // Here we could
405 // a) check the type
406 //fprintf(stderr,"found #%d arg %s\n",nargs,arg.Data());
407
408 // We register the arg for later usage
409 argArr.Add(new TObjString(arg));
410 nargs++;
411
412 prevComma = i+1;
413 };
414 }
415
416 if (nargs>999) {
417 err = 7;
418 return false;
419 }
420
421 // Now we need to lookup the function and check its arguments.
422 TClass *ns = (spaceName.Length()) ? TClass::GetClass(spaceName) : 0;
423 ClassInfo_t *cinfo = 0;
424 if (ns) {
425 cinfo = ns->GetClassInfo();
426 } else {
427 cinfo = gInterpreter->ClassInfo_Factory();
428 }
429
430 // ROOT does yet have a complete TType class, but TCling does,
431 // so let's use that for now.
432 static TypeInfo_t *const doubletype { gInterpreter->TypeInfo_Factory("double") };
433
434 std::vector<TypeInfo_t*> proto(nargs,doubletype);
435
436 CallFunc_t *callfunc = gInterpreter->CallFunc_Factory();
437 Long_t func_offset;
438 gInterpreter->CallFunc_SetFuncProto(callfunc,cinfo,functionName,proto,false,&func_offset,ROOT::kConversionMatch);
439
440 TMethodCall *method = new TMethodCall(ns,callfunc,func_offset);
441
442 if (!ns) gInterpreter->ClassInfo_Delete(cinfo);
443 gInterpreter->CallFunc_Delete(callfunc);
444
445 if (method->IsValid()) {
446 if (method->ReturnType() == TMethodCall::kOther) {
447 /*
448 Error("Compile",
449 "TFormula can only call interpreted and compiled function that returns a numerical type %s returns a %s\n",
450 method->GetMethodName(), method->GetMethod()->GetReturnTypeName());
451 */
452 err=29;
453
454 } else {
455
456 // Analyze the arguments
457 TIter next(&argArr);
458 TObjString *objstr;
459 while ( (objstr=(TObjString*)next()) ) {
460 Analyze(objstr->String(),err,offset);
461 }
462
463 fFunctions.Add(method);
464 fExpr[fNoper] = method->GetMethod()->GetPrototype();
465 SetAction(fNoper, kFunctionCall, fFunctions.GetLast()*1000 + nargs);
466 fNoper++;
467 return true;
468 }
469 }
470
471 delete method;
472 //
473 // MI change - extended space of functions
474 // not forward compatible change
475 //
476 TString cbase(chaine);
477 Int_t args_paran = cbase.First("(");
478 if (args_paran>0){
479 cbase[args_paran]=0;
480 }
481
482 ROOT::v5::TFormulaPrimitive *prim = ROOT::v5::TFormulaPrimitive::FindFormula(cbase, args_paran>0 ? cbase.Data() + args_paran + 1 : (const char*)0);
483 if (prim && (!IsA()->GetBaseClass("TTreeFormula"))) {
484 // TO BE DONE ALSO IN TTREFORMULA - temporary fix MI
485 // Analyze the arguments
486 TIter next(&argArr);
487 TObjString *objstr;
488 while ( (objstr=(TObjString*)next()) ) {
489 Analyze(objstr->String(),err,offset); if (err) return kFALSE;
490 }
491 if (nargs!=prim->fNArguments) {
492 Error("Compile", "%s requires %d arguments",
493 prim->GetName(), prim->fNArguments);
494 return kFALSE;
495 }
496 fExpr[fNoper] = prim->GetName();
497 if (prim->fType==10){
498 SetAction(fNoper, kFD1);
499 }
500 if (prim->fType==110){
501 SetAction(fNoper, kFD2);
502 }
503 if (prim->fType==1110){
504 SetAction(fNoper, kFD3);
505 }
506 if (prim->fType==-1){
507 SetAction(fNoper, kFDM);
508 if (fNpar<prim->fNParameters) fNpar+=prim->fNParameters;
509 }
510
511 fNoper++;
512 return kTRUE;
513 }
514
515 return kFALSE;
516}
517
518////////////////////////////////////////////////////////////////////////////////
519/// Analyze a sub-expression in one formula.
520///
521/// Expressions in one formula are recursively analyzed.
522/// Result of analysis is stored in the object tables.
523///
524/// ### Table of function codes and errors
525///
526/// ~~~ {.cpp}
527/// * functions :
528///
529/// + 1 pow 20
530/// - 2 sq 21
531/// * 3 sqrt 22
532/// / 4 strstr 23
533/// % 5 min 24
534/// max 25
535/// log 30
536/// cos 10 exp 31
537/// sin 11 log10 32
538/// tan 12
539/// acos 13 abs 41
540/// asin 14 sign 42
541/// atan 15 int 43
542/// atan2 16
543/// fmod 17 rndm 50
544///
545/// cosh 70 acosh 73
546/// sinh 71 asinh 74
547/// tanh 72 atanh 75
548///
549/// expo 100 gaus 110 gausn (see note below)
550/// expo(0) 100 0 gaus(0) 110 0 gausn(0)
551/// expo(1) 100 1 gaus(1) 110 1 gausn(1)
552/// xexpo 100 x xgaus 110 x xgausn
553/// yexpo 101 x ygaus 111 x ygausn
554/// zexpo 102 x zgaus 112 x zgausn
555/// xyexpo 105 x xygaus 115 x xygausn
556/// yexpo(5) 102 5 ygaus(5) 111 5 ygausn(5)
557/// xyexpo(2) 105 2 xygaus(2) 115 2 xygausn(2)
558///
559/// landau 120 x landaun (see note below)
560/// landau(0) 120 0 landaun(0)
561/// landau(1) 120 1 landaun(1)
562/// xlandau 120 x xlandaun
563/// ylandau 121 x ylandaun
564/// zlandau 122 x zlandaun
565/// xylandau 125 x xylandaun
566/// ylandau(5) 121 5 ylandaun(5)
567/// xylandau(2) 125 2 xylandaun(2)
568///
569/// pol0 130 x pol1 130 1xx
570/// pol0(0) 130 0 pol1(0) 130 100
571/// pol0(1) 130 1 pol1(1) 130 101
572/// xpol0 130 x xpol1 130 101
573/// ypol0 131 x ypol1 131 101
574/// zpol0 132 x zpol1 132 1xx
575/// ypol0(5) 131 5 ypol1(5) 131 105
576///
577/// pi 40
578///
579/// && 60 < 64
580/// || 61 > 65
581/// == 62 <= 66
582/// != 63 => 67
583/// ! 68
584/// ==(string) 76 & 78
585/// !=(string) 77 | 79
586/// <<(shift) 80 >>(shift) 81
587/// ? : 82
588///
589/// * constants (kConstants) :
590///
591/// c0 141 1 c1 141 2 etc..
592///
593/// * strings (kStringConst):
594///
595/// sX 143 x
596///
597/// * variables (kFormulaVar) :
598///
599/// x 144 0 y 144 1 z 144 2 t 144 3
600///
601/// * parameters :
602///
603/// [1] 140 1
604/// [2] 140 2
605/// etc.
606/// ~~~
607///
608/// ### Special cases for normalized gaussian or landau distributions
609///
610/// the expression "gaus" is a substitute for
611///
612/// [0]*exp(-0.5*((x-[1])/[2])**2)
613///
614/// to obtain a standard normalized gaussian, use "gausn" instead of "gaus"
615/// the expression "gausn" is a substitute for
616///
617/// [0]*exp(-0.5*((x-[1])/[2])**2)/(sqrt(2*pi)*[2]))
618///
619/// WARNING: gaus and gausn are mutually exclusive in the same expression.
620///
621/// In the same way the expression "landau" is a substitute for
622///
623/// [0]*TMath::Landau(x,[1],[2],kFALSE)
624///
625/// to obtain a standard normalized landau, use "landaun" instead of "landau"
626/// the expression "landaun" is a substitute for
627///
628/// [0]*TMath::Landau(x,[1],[2],kTRUE)
629///
630/// WARNING: landau and landaun are mutually exclusive in the same expression.
631///
632/// ### Boolean optimization (kBoolOptmize) :
633///
634/// Those pseudo operation are used to implement lazy evaluation of
635/// && and ||. When the left hand of the expression if false
636/// (respectively true), the evaluation of the right is entirely skipped
637/// (since it would not change the value of the expression).
638///
639/// && 142 11 (one operation on right) 142 21 (2 operations on right)
640/// || 142 12 (one operation on right) 142 22 (2 operations on right)
641///
642/// * functions calls (kFunctionCall) :
643///
644/// f0 145 0 f1 145 1 etc..
645///
646/// ### Errors :
647///
648/// 1 : Division By Zero
649/// 2 : Invalid Floating Point Operation
650/// 4 : Empty String
651/// 5 : invalid syntax
652/// 6 : Too many operators
653/// 7 : Too many parameters
654/// 10 : z specified but not x and y
655/// 11 : z and y specified but not x
656/// 12 : y specified but not x
657/// 13 : z and x specified but not y
658/// 20 : non integer value for parameter number
659/// 21 : atan2 requires two arguments
660/// 22 : pow requires two arguments
661/// 23 : degree of polynomial not specified
662/// 24 : Degree of polynomial must be positive
663/// 25 : Degree of polynomial must be less than 20
664/// 26 : Unknown name
665/// 27 : Too many constants in expression
666/// 28 : strstr requires two arguments
667/// 29 : interpreted or compiled function have to return a numerical type
668/// 30 : Bad numerical expression
669/// 31 : Part of the variable exist but some of it is not accessible or useable
670/// 40 : '(' is expected
671/// 41 : ')' is expected
672/// 42 : '[' is expected
673/// 43 : ']' is expected
674///
675/// \image html TFormula_analyze.png
676///
677/// ### Special functions
678///
679/// By default, the formula is assigned fNumber=0. However, the following
680/// formula built with simple functions are assigned fNumber:
681///
682/// "gaus" 100 (or gausn)
683/// "xygaus" 110
684/// "expo" 200
685/// "polN" 300+N
686/// "landau" 400
687/// "xylandau" 410
688///
689/// Note that expressions like gaus(0), expo(1) will force fNumber=0
690///
691/// ### Warning when deriving a class from TFormula
692///
693/// If you overload this member function, you also HAVE TO
694/// never call the constructor:
695///
696/// ~~~ {.cpp}
697/// TFormula::TFormula(const char *name,const char *expression)
698/// ~~~
699///
700/// and write your own constructor
701///
702/// ~~~ {.cpp}
703/// MyClass::MyClass(const char *name,const char *expression) : TFormula()
704/// ~~~
705///
706/// which has to call the TFormula default constructor and whose implementation
707/// should be similar to the implementation of the normal TFormula constructor
708///
709/// This is necessary because the normal TFormula constructor call indirectly
710/// the virtual member functions Analyze, DefaultString, DefaultValue
711/// and DefaultVariable.
712
713void TFormula::Analyze(const char *schain, Int_t &err, Int_t offset)
714{
715
716
717 Int_t valeur,find,n,i,j,k,lchain,nomb,virgule,inter,nest;
718 valeur=find=n=i=j=k=lchain=nomb=virgule=inter=nest = 0;
719 Int_t compt,compt2,compt3,compt4;
720 Bool_t inString;
721 Double_t vafConst;
722 ULong_t vafConst2;
723 Bool_t parenthese;
724 TString s,chaine_error,chaine1ST;
725 TString s1,s2,s3,ctemp;
726
727 TString chaine = schain;
728 const TFormula *oldformula;
729 Int_t modulo,plus,puiss10,puiss10bis,moins,multi,divi,puiss,et,ou,petit,grand,egal,diff,peteg,grdeg,etx,oux,rshift,lshift,tercond,terelse;
730 char t;
731 TString slash("/"), escapedSlash("\\/");
732 Int_t inter2 = 0;
733 SetNumber(0);
734 Int_t actionCode,actionParam;
735 Int_t err_hint = 0;
736
737 // Verify correct matching of parenthesis and remove unnecessary parenthesis.
738 lchain = chaine.Length();
739 //if (chaine(lchain-2,2) == "^2") chaine = "sq(" + chaine(0,lchain-2) + ")";
740 parenthese = kTRUE;
741 lchain = chaine.Length();
742 while (parenthese && lchain>0 && err==0){
743 compt = 0;
744 compt2 = 0;
745 inString = false;
746 lchain = chaine.Length();
747 if (lchain==0) err=4;
748 else {
749 for (i=1; i<=lchain; ++i) {
750 if (chaine(i-1,1) == "\"") inString = !inString;
751 if (!inString) {
752 if (chaine(i-1,1) == "[") compt2++;
753 if (chaine(i-1,1) == "]") compt2--;
754 if (chaine(i-1,1) == "(") compt++;
755 if (chaine(i-1,1) == ")") compt--;
756 }
757 if (compt < 0) err = 40; // more open parentheses than close parentheses
758 if (compt2< 0) err = 42; // more ] than [
759 if (compt==0 && (i!=lchain || lchain==1)) parenthese = kFALSE;
760 // if (lchain<3 && chaine(0,1)!="(" && chaine(lchain-1,1)!=")") parenthese = kFALSE;
761 }
762 if (compt > 0) err = 41; // more ( than )
763 if (compt2> 0) err = 43; // more [ than ]
764 if (parenthese) chaine = chaine(1,lchain-2);
765 }
766 } // while parenthesis
767
768 if (lchain==0) err=4; // empty string
769 modulo=plus=moins=multi=divi=puiss=et=ou=petit=grand=egal=diff=peteg=grdeg=etx=oux=rshift=lshift=tercond=terelse=0;
770
771 // Look for simple operators
772
773 if (err==0) {
774 compt = compt2 = compt3 = compt4 = 0;puiss10=0;puiss10bis = 0;
775 inString = false;
776 j = lchain;
777 Bool_t isdecimal = 1; // indicates whether the left part is decimal.
778
779 for (i=1;i<=lchain; i++) {
780
781 puiss10=puiss10bis=0;
782 if (i>2) {
783 t = chaine[i-3];
784 isdecimal = isdecimal && (strchr("0123456789.",t)!=0);
785 if (isdecimal) {
786 if ( chaine[i-2] == 'e' || chaine[i-2] == 'E' ) puiss10 = 1;
787 } else if ( strchr("+-/[]()&|><=!*/%^\\",t) ) {
788 isdecimal = 1; // reset after delimiter
789 }
790 }
791 if (j>2) {
792 if (chaine[j-2] == 'e' || chaine[j-2] == 'E') {
793 Bool_t isrightdecimal = 1;
794
795 for(k=j-3; k>=0 && isrightdecimal; --k) {
796 t = chaine[k];
797 isrightdecimal = isrightdecimal && (strchr("0123456789.",t)!=0);
798 if (!isrightdecimal) {
799 if (strchr("+-/[]()&|><=!*/%^\\",t)!=0) {
800 puiss10bis = 1;
801 }
802 }
803 }
804 if (k<0 && isrightdecimal) puiss10bis = 1;
805 }
806 }
807 if (puiss10 && (i<=lchain)) {
808 t = chaine[i];
809 puiss10 = (strchr("0123456789.",t)!=0);
810 }
811 if (puiss10bis && (j<=lchain)) {
812 t = chaine[j];
813 puiss10bis = (strchr("0123456789.",t)!=0);
814 }
815
816 if (chaine(i-1,1) == "\"") inString = !inString;
817 if (inString) continue;
818 if (chaine(i-1,1) == "[") compt2++;
819 if (chaine(i-1,1) == "]") compt2--;
820 if (chaine(i-1,1) == "(") compt++;
821 if (chaine(i-1,1) == ")") compt--;
822 if (chaine(j-1,1) == "[") compt3++;
823 if (chaine(j-1,1) == "]") compt3--;
824 if (chaine(j-1,1) == "(") compt4++;
825 if (chaine(j-1,1) == ")") compt4--;
826 if (chaine(i-1,2)=="&&" && !inString && compt==0 && compt2==0 && et==0) {et=i;puiss=0;}
827 if (chaine(i-1,2)=="||" && compt==0 && compt2==0 && ou==0) {puiss10=0; ou=i;}
828 if (chaine(i-1,1)=="&" && compt==0 && compt2==0 && etx==0) {etx=i;puiss=0;}
829 if (chaine(i-1,1)=="|" && compt==0 && compt2==0 && oux==0) {puiss10=0; oux=i;}
830 if (chaine(i-1,2)==">>" && compt==0 && compt2==0 && rshift==0) {puiss10=0; rshift=i;}
831 if (chaine(i-1,1)==">" && compt==0 && compt2==0 && rshift==0 && grand==0)
832 {puiss10=0; grand=i;}
833 if (chaine(i-1,2)=="<<" && compt==0 && compt2==0 && lshift==0) {puiss10=0; lshift=i;}
834 if (chaine(i-1,1)=="<" && compt==0 && compt2==0 && lshift==0 && petit==0)
835 {puiss10=0; petit=i;
836 // Check whether or not we have a template names! (actually this can
837 // only happen in TTreeFormula.
838 for(int ip = i,depth=0; ip < lchain; ++ip) {
839 char c = chaine(ip);
840 // The characters allowed in the template parameter are alpha-numerical characters,
841 // underscores, comma, <, > and scope operator.
842 if (isalnum(c) || c=='_' || c==',') continue;
843 if (c==':' && chaine(ip+1)==':') { ++ip; continue; }
844 if (c=='<') { ++depth; continue; }
845 if (c=='>') {
846 if (depth) { --depth; continue; }
847 else {
848 // We reach the end of the template parameter.
849 petit = 0;
850 i = ip+1;
851 break;
852 }
853 }
854 // Character not authorized within a template parameter
855 break;
856 }
857 if (petit==0) {
858 // We found a template parameter and modified i
859 continue; // the for(int i ,...)
860 }
861 }
862 if ((chaine(i-1,2)=="<=" || chaine(i-1,2)=="=<") && compt==0 && compt2==0
863 && peteg==0) {peteg=i; puiss10=0; petit=0;}
864 if ((chaine(i-1,2)=="=>" || chaine(i-1,2)==">=") && compt==0 && compt2==0
865 && grdeg==0) {puiss10=0; grdeg=i; grand=0;}
866 if (chaine(i-1,2) == "==" && compt == 0 && compt2 == 0 && egal == 0) {puiss10=0; egal=i;}
867 if (chaine(i-1,2) == "!=" && compt == 0 && compt2 == 0 && diff == 0) {puiss10=0; diff=i;}
868 if (i>1 && chaine(i-1,1) == "+" && compt == 0 && compt2 == 0 && puiss10==0) plus=i;
869 if (chaine(j-1,1) == "-" && chaine(j-2,1) != "*" && chaine(j-2,1) != "/"
870 && chaine(j-2,1)!="^" && compt3==0 && compt4==0 && moins==0 && puiss10bis==0) moins=j;
871 if (chaine(i-1,1)=="%" && compt==0 && compt2==0 && modulo==0) {puiss10=0; modulo=i;}
872 if (chaine(i-1,1)=="*" && compt==0 && compt2==0 && multi==0) {puiss10=0; multi=i;}
873 if (chaine(j-1,1)=="/" && chaine(j-2,1)!="\\"
874 && compt4==0 && compt3==0 && divi==0)
875 {
876 puiss10=0; divi=j;
877 }
878 if (chaine(j-1)=='^' && compt4==0 && compt3==0 && puiss==0) {puiss10=0; puiss=j;}
879 if (chaine(i-1)=='?' && compt == 0 && compt2 == 0 && tercond == 0) {puiss10=0; tercond=i;}
880 if (chaine(i-1)==':' && tercond && compt == 0 && compt2 == 0 && terelse == 0) {
881 if (i>2 && chaine(i-2)!=':' && chaine(i)!=':') {
882 puiss10=0; terelse=i;
883 }
884 }
885
886 j--;
887 }
888
889 // If operator found, analyze left and right part of the statement
890
891 actionParam = 0;
892 if (tercond && terelse) {
893 if (tercond == 1 || terelse == lchain || tercond == (terelse-1) ) {
894 err = 5;
895 chaine_error = "?:";
896 } else {
897 // Condition
898 ctemp = chaine(0,tercond-1);
899 Analyze(ctemp.Data(),err,offset); if (err) return;
900
901 fExpr[fNoper] = "?: condition jump";
902 actionCode = kJumpIf;
903 actionParam = 0;
904 SetAction(fNoper,actionCode, actionParam);
905 Int_t optloc = fNoper++;
906
907 // Expression executed if condition is true.
908 ctemp = chaine(tercond,terelse-tercond-1);
909 Analyze(ctemp.Data(),err,offset); if (err) return;
910 actionParam = fNoper; // We want to skip the next instruction ('else jump'), so we set the param to the current cursor and the next instruction will be skip by the ++i in the eval loop
911 SetAction(optloc, actionCode, actionParam);
912
913 fExpr[fNoper] = "?: else jump";
914 actionCode = kJump;
915 actionParam = 0;
916 // Set jump target.
917 SetAction(fNoper,actionCode, actionParam);
918 optloc = fNoper++;
919
920 // Expression executed if condition is false.
921 ctemp = chaine(terelse,lchain-terelse);
922 Analyze(ctemp.Data(),err,offset); if (err) return;
923 // Set jump target.
924 actionParam = fNoper - 1; // We need to not skip the next instruction, so we compensate for the ++i in the eval loop
925 SetAction(optloc, actionCode, actionParam);
926
927 if (IsString(optloc-1) != IsString(fNoper-1)) {
928 err = 45;
929 chaine_error = "?:";
930 }
931 }
932 } else if (ou != 0) { //check for ||
933 if (ou==1 || ou==lchain-1) {
934 err=5;
935 chaine_error="||";
936 }
937 else {
938 ctemp = chaine(0,ou-1);
939 Analyze(ctemp.Data(),err,offset); if (err) return;
940
941 fExpr[fNoper] = "|| checkpoint";
942 actionCode = kBoolOptimize;
943 actionParam = 2;
944 SetAction(fNoper,actionCode, actionParam);
945 Int_t optloc = fNoper++;
946
947 ctemp = chaine(ou+1,lchain-ou-1);
948 Analyze(ctemp.Data(),err,offset); if (err) return;
949 fExpr[fNoper] = "||";
950 actionCode = kOr;
951 SetAction(fNoper,actionCode, 0);
952
953 SetAction( optloc, GetAction(optloc), GetActionParam(optloc) + (fNoper-optloc) * 10);
954 fNoper++;
955 if (!CheckOperands(optloc-1,fNoper-1,err)) return;
956 }
957 } else if (et!=0) {
958 if (et==1 || et==lchain-1) {
959 err=5;
960 chaine_error="&&";
961 }
962 else {
963 ctemp = chaine(0,et-1);
964 Analyze(ctemp.Data(),err,offset); if (err) return;
965
966 fExpr[fNoper] = "&& checkpoint";
967 actionCode = kBoolOptimize;
968 actionParam = 1;
969 SetAction(fNoper,actionCode,actionParam);
970
971 Int_t optloc = fNoper++;
972
973 ctemp = chaine(et+1,lchain-et-1);
974 Analyze(ctemp.Data(),err,offset); if (err) return;
975 fExpr[fNoper] = "&&";
976 actionCode = kAnd;
977 SetAction(fNoper,actionCode,0);
978
979 SetAction(optloc, GetAction(optloc), GetActionParam(optloc) + (fNoper-optloc) * 10);
980 fNoper++;
981 if (!CheckOperands(optloc-1,fNoper-1,err)) return;
982 }
983 } else if (oux!=0) {
984 if (oux==1 || oux==lchain) {
985 err=5;
986 chaine_error="|";
987 }
988 else {
989 ctemp = chaine(0,oux-1);
990 Analyze(ctemp.Data(),err,offset); if (err) return;
991 UInt_t leftopr = fNoper-1;
992 ctemp = chaine(oux,lchain-oux);
993 Analyze(ctemp.Data(),err,offset); if (err) return;
994 fExpr[fNoper] = "|";
995 actionCode = kBitOr;
996 SetAction(fNoper,actionCode,actionParam);
997 fNoper++;
998 if (!CheckOperands(leftopr,fNoper-1,err)) return;
999 }
1000 } else if (etx!=0) {
1001 if (etx==1 || etx==lchain) {
1002 err=5;
1003 chaine_error="&";
1004 }
1005 else {
1006 ctemp = chaine(0,etx-1);
1007 Analyze(ctemp.Data(),err,offset); if (err) return;
1008 UInt_t leftopr = fNoper-1;
1009 ctemp = chaine(etx,lchain-etx);
1010 Analyze(ctemp.Data(),err,offset); if (err) return;
1011 fExpr[fNoper] = "&";
1012 actionCode = kBitAnd;
1013 SetAction(fNoper,actionCode,actionParam);
1014 fNoper++;
1015 if (!CheckOperands(leftopr,fNoper-1,err)) return;
1016 }
1017 } else if (petit != 0) {
1018 if (petit==1 || petit==lchain) {
1019 err=5;
1020 chaine_error="<";
1021 }
1022 else {
1023 ctemp = chaine(0,petit-1);
1024 Analyze(ctemp.Data(),err,offset); if (err) return;
1025 UInt_t leftopr = fNoper-1;
1026 ctemp = chaine(petit,lchain-petit);
1027 Analyze(ctemp.Data(),err,offset); if (err) return;
1028 fExpr[fNoper] = "<";
1029 actionCode = kLess;
1030 SetAction(fNoper,actionCode,actionParam);
1031 fNoper++;
1032 if (!CheckOperands(leftopr,fNoper-1,err)) return;
1033 }
1034 } else if (grand != 0) {
1035 if (grand==1 || grand==lchain) {
1036 err=5;
1037 chaine_error=">";
1038 }
1039 else {
1040 ctemp = chaine(0,grand-1);
1041 Analyze(ctemp.Data(),err,offset); if (err) return;
1042 UInt_t leftopr = fNoper-1;
1043 ctemp = chaine(grand,lchain-grand);
1044 Analyze(ctemp.Data(),err,offset); if (err) return;
1045 fExpr[fNoper] = ">";
1046 actionCode = kGreater;
1047 SetAction(fNoper,actionCode,actionParam);
1048 fNoper++;
1049 if (!CheckOperands(leftopr,fNoper-1,err)) return;
1050 }
1051 } else if (peteg != 0) {
1052 if (peteg==1 || peteg==lchain-1) {
1053 err=5;
1054 chaine_error="<=";
1055 }
1056 else {
1057 ctemp = chaine(0,peteg-1);
1058 Analyze(ctemp.Data(),err,offset); if (err) return;
1059 ctemp = chaine(peteg+1,lchain-peteg-1);
1060 UInt_t leftopr = fNoper-1;
1061 Analyze(ctemp.Data(),err,offset); if (err) return;
1062 fExpr[fNoper] = "<=";
1063 actionCode = kLessThan;
1064 SetAction(fNoper,actionCode,actionParam);
1065 fNoper++;
1066 if (!CheckOperands(leftopr,fNoper-1,err)) return;
1067 }
1068 } else if (grdeg != 0) {
1069 if (grdeg==1 || grdeg==lchain-1) {
1070 err=5;
1071 chaine_error="=>";
1072 }
1073 else {
1074 ctemp = chaine(0,grdeg-1);
1075 Analyze(ctemp.Data(),err,offset); if (err) return;
1076 UInt_t leftopr = fNoper-1;
1077 ctemp = chaine(grdeg+1,lchain-grdeg-1);
1078 Analyze(ctemp.Data(),err,offset); if (err) return;
1079 fExpr[fNoper] = ">=";
1080 actionCode = kGreaterThan;
1081 SetAction(fNoper,actionCode,actionParam);
1082 fNoper++;
1083 if (!CheckOperands(leftopr,fNoper-1,err)) return;
1084 }
1085 } else if (egal != 0) {
1086 if (egal==1 || egal==lchain-1) {
1087 err=5;
1088 chaine_error="==";
1089 }
1090 else {
1091 ctemp = chaine(0,egal-1);
1092 Analyze(ctemp.Data(),err,offset); if (err) return;
1093 Int_t optloc = fNoper-1;
1094
1095 ctemp = chaine(egal+1,lchain-egal-1);
1096 Analyze(ctemp.Data(),err,offset); if (err) return;
1097 fExpr[fNoper] = "==";
1098 actionCode = kEqual;
1099
1100 Bool_t isstring = IsString(fNoper-1);
1101 if (IsString(optloc) != isstring) {
1102 err = 45;
1103 chaine_error = "==";
1104 } else if (isstring) {
1105 actionCode = kStringEqual;
1106 }
1107 SetAction(fNoper,actionCode,actionParam);
1108 fNoper++;
1109 }
1110 } else if (diff != 0) {
1111 if (diff==1 || diff==lchain-1) {
1112 err=5;
1113 chaine_error = "!=";
1114 }
1115 else {
1116 ctemp = chaine(0,diff-1);
1117 Analyze(ctemp.Data(),err,offset); if (err) return;
1118 Int_t optloc = fNoper-1;
1119
1120 ctemp = chaine(diff+1,lchain-diff-1);
1121 Analyze(ctemp.Data(),err,offset); if (err) return;
1122 fExpr[fNoper] = "!=";
1123 actionCode = kNotEqual;
1124
1125 Bool_t isstring = IsString(fNoper-1);
1126 if (IsString(optloc) != isstring) {
1127 err = 45;
1128 chaine_error = "!=";
1129 } else if (isstring) {
1130 actionCode = kStringNotEqual;
1131 }
1132 SetAction(fNoper,actionCode,actionParam);
1133 fNoper++;
1134 }
1135 } else if (plus != 0) {
1136 if (plus==lchain) {
1137 err=5;
1138 chaine_error = "+";
1139 }
1140 else {
1141 ctemp = chaine(0,plus-1);
1142 Analyze(ctemp.Data(),err,offset); if (err) return;
1143 UInt_t leftopr = fNoper-1;
1144 ctemp = chaine(plus,lchain-plus);
1145 Analyze(ctemp.Data(),err,offset); if (err) return;
1146 fExpr[fNoper] = "+";
1147 actionCode = kAdd;
1148 SetAction(fNoper,actionCode,actionParam);
1149 fNoper++;
1150 if (!CheckOperands(leftopr,fNoper-1,err)) return;
1151 }
1152 } else {
1153 if (moins != 0) {
1154 if (moins == 1) {
1155 ctemp = chaine(moins,lchain-moins);
1156 Analyze(ctemp.Data(),err,offset); if (err) return;
1157 fExpr[fNoper] = "-";
1158 actionCode = kSignInv;
1159 SetAction(fNoper,actionCode,actionParam);
1160 ++fNoper;
1161 if (!CheckOperands(fNoper-1,err)) return;
1162 } else {
1163 if (moins == lchain) {
1164 err=5;
1165 chaine_error = "-";
1166 } else {
1167 ctemp = chaine(0,moins-1);
1168 Analyze(ctemp.Data(),err,offset); if (err) return;
1169 UInt_t leftopr = fNoper-1;
1170 ctemp = chaine(moins,lchain-moins);
1171 Analyze(ctemp.Data(),err,offset); if (err) return;
1172 fExpr[fNoper] = "-";
1173 actionCode = kSubstract;
1174 SetAction(fNoper,actionCode,actionParam);
1175 fNoper++;
1176 if (!CheckOperands(leftopr,fNoper-1,err)) return;
1177 }
1178 }
1179 } else if (modulo != 0) {
1180 if (modulo == 1 || modulo == lchain) {
1181 err=5;
1182 chaine_error="%";
1183 } else {
1184 ctemp = chaine(0,modulo-1);
1185 Analyze(ctemp.Data(),err,offset); if (err) return;
1186 UInt_t leftopr = fNoper-1;
1187 ctemp = chaine(modulo,lchain-modulo);
1188 Analyze(ctemp.Data(),err,offset); if (err) return;
1189 fExpr[fNoper] = "%";
1190 actionCode = kModulo;
1191 SetAction(fNoper,actionCode,actionParam);
1192 fNoper++;
1193 if (!CheckOperands(leftopr,fNoper-1,err)) return;
1194 }
1195 } else if (rshift != 0) {
1196 if (rshift == 1 || rshift == lchain) {
1197 err=5;
1198 chaine_error=">>";
1199 } else {
1200 ctemp = chaine(0,rshift-1);
1201 Analyze(ctemp.Data(),err,offset); if (err) return;
1202 UInt_t leftopr = fNoper-1;
1203 ctemp = chaine(rshift+1,lchain-rshift-1);
1204 Analyze(ctemp.Data(),err,offset); if (err) return;
1205 fExpr[fNoper] = ">>";
1206 actionCode = kRightShift;
1207 SetAction(fNoper,actionCode,actionParam);
1208 fNoper++;
1209 if (!CheckOperands(leftopr,fNoper-1,err)) return;
1210 }
1211 } else if (lshift != 0) {
1212 if (lshift == 1 || lshift == lchain) {
1213 err=5;
1214 chaine_error=">>";
1215 } else {
1216 ctemp = chaine(0,lshift-1);
1217 Analyze(ctemp.Data(),err,offset); if (err) return;
1218 UInt_t leftopr = fNoper-1;
1219 ctemp = chaine(lshift+1,lchain-lshift-1);
1220 Analyze(ctemp.Data(),err,offset); if (err) return;
1221 fExpr[fNoper] = ">>";
1222 actionCode = kLeftShift;
1223 SetAction(fNoper,actionCode,actionParam);
1224 fNoper++;
1225 if (!CheckOperands(leftopr,fNoper-1,err)) return;
1226 }
1227 } else {
1228 if (multi != 0) {
1229 if (multi == 1 || multi == lchain) {
1230 err=5;
1231 chaine_error="*";
1232 }
1233 else {
1234 ctemp = chaine(0,multi-1);
1235 Analyze(ctemp.Data(),err,offset); if (err) return;
1236 UInt_t leftopr = fNoper-1;
1237 ctemp = chaine(multi,lchain-multi);
1238 Analyze(ctemp.Data(),err,offset); if (err) return;
1239 fExpr[fNoper] = "*";
1240 actionCode = kMultiply;
1241 SetAction(fNoper,actionCode,actionParam);
1242 fNoper++;
1243 if (!CheckOperands(leftopr,fNoper-1,err)) return;
1244 }
1245 } else {
1246 if (divi != 0) {
1247 if (divi == 1 || divi == lchain) {
1248 err=5;
1249 chaine_error = "/";
1250 }
1251 else {
1252 ctemp = chaine(0,divi-1);
1253 Analyze(ctemp.Data(),err,offset); if (err) return;
1254 UInt_t leftopr = fNoper-1;
1255 ctemp = chaine(divi,lchain-divi);
1256 Analyze(ctemp.Data(),err,offset); if (err) return;
1257 fExpr[fNoper] = "/";
1258 actionCode = kDivide;
1259 SetAction(fNoper,actionCode,actionParam);
1260 fNoper++;
1261 if (!CheckOperands(leftopr,fNoper-1,err)) return;
1262 }
1263 } else {
1264 if (puiss != 0) {
1265 if (puiss == 1 || puiss == lchain) {
1266 err = 5;
1267 chaine_error = "**";
1268 }
1269 else {
1270 if (chaine(lchain-2,2) == "^2") {
1271 ctemp = "sq(" + chaine(0,lchain-2) + ")";
1272 Analyze(ctemp.Data(),err,offset); if (err) return;
1273 } else {
1274 ctemp = chaine(0,puiss-1);
1275 Analyze(ctemp.Data(),err,offset); if (err) return;
1276 UInt_t leftopr = fNoper-1;
1277 ctemp = chaine(puiss,lchain-puiss);
1278 Analyze(ctemp.Data(),err,offset); if (err) return;
1279 fExpr[fNoper] = "^";
1280 actionCode = kpow;
1281 SetAction(fNoper,actionCode,actionParam);
1282 fNoper++;
1283 if (!CheckOperands(leftopr,fNoper-1,err)) return;
1284 }
1285 }
1286 } else {
1287
1288 find=0;
1289
1290 // Check for a numerical expression
1291 {
1292 Bool_t hasDot = kFALSE;
1293 Bool_t isHexa = kFALSE;
1294 Bool_t hasExpo= kFALSE;
1295 if ((chaine(0,2)=="0x")||(chaine(0,2)=="0X")) isHexa=kTRUE;
1296 for (j=0; j<chaine.Length() && err==0; j++) {
1297 t=chaine[j];
1298 if (!isHexa) {
1299 if (j>0 && (chaine(j,1)=="e" || chaine(j,2)=="e+" || chaine(j,2)=="e-" || chaine(j,1)=="E" || chaine(j,2)=="E+" || chaine(j,2)=="E-")) {
1300 if (hasExpo) {
1301 err=26;
1302 chaine_error=chaine;
1303 }
1304 hasExpo = kTRUE;
1305 // The previous implementation allowed a '.' in the exponent.
1306 // That information was ignored (by sscanf), we now make it an error
1307 // hasDot = kFALSE;
1308 hasDot = kTRUE; // forbid any additional '.'
1309 if (chaine(j,2)=="e+" || chaine(j,2)=="e-" || chaine(j,2)=="E+" || chaine(j,2)=="E-") j++;
1310 }
1311 else {
1312 if (chaine(j,1) == "." && !hasDot) hasDot = kTRUE; // accept only one '.' in the number
1313 else {
1314 // The previous implementation was allowing ANYTHING after the '.' and thus
1315 // made legal code like '2.3 and fpx' and was just silently ignoring the
1316 // 'and fpx'.
1317 if (!strchr("0123456789",t) && (chaine(j,1)!="+" || j!=0)) {
1318 err = 30;
1319 chaine_error=chaine;
1320 }
1321 }
1322 }
1323 }
1324 else {
1325 if (!strchr("0123456789abcdefABCDEF",t) && (j>1)) {
1326 err = 30;
1327 chaine_error=chaine;
1328 }
1329 }
1330 }
1331 if (fNconst >= gMAXCONST) err = 27;
1332 if (!err) {
1333 if (!isHexa) {if (sscanf((const char*)chaine,"%lg",&vafConst) > 0) err = 0; else err =1;}
1334 else {if (sscanf((const char*)chaine,"%lx",&vafConst2) > 0) err = 0; else err=1;
1335 vafConst = (Double_t) vafConst2;}
1336 fExpr[fNoper] = chaine;
1337 k = -1;
1338 for (j=0;j<fNconst;j++) {
1339 if (vafConst == fConst[j] ) k= j;
1340 }
1341 if ( k < 0) { k = fNconst; fNconst++; fConst[k] = vafConst; }
1342 actionCode = kConstant;
1343 actionParam = k;
1344 SetAction(fNoper,actionCode,actionParam);
1345 fNoper++;
1346 }
1347 if (err==30) err=0;
1348 else find = kTRUE;
1349 }
1350
1351 // Look for an already defined expression
1352
1353 if (find==0) {
1354 {
1355 R__LOCKGUARD(gROOTMutex);
1356 oldformula = (const TFormula*)gROOT->GetListOfFunctions()->FindObject((const char*)chaine);
1357 }
1358 if (oldformula && strcmp(schain,oldformula->GetTitle())) {
1359 Int_t nprior = fNpar;
1360 Analyze(oldformula->GetExpFormula(),err,fNpar);
1361 // if the oldformula was using a normalized function (gausn or landaun) set also in this one
1362 if (oldformula->IsNormalized()) SetBit(kNormalized);
1363 if (err) return; // changes fNpar
1364 fNpar = nprior;
1365 find=1;
1366 if (!err) {
1367 Int_t npold = oldformula->GetNpar();
1368 fNpar += npold;
1369 for (Int_t ipar=0;ipar<npold;ipar++) {
1370 fParams[ipar+fNpar-npold] = oldformula->GetParameter(ipar);
1371 }
1372 }
1373 }
1374 }
1375 if (find == 0) {
1376
1377 // Check if chaine is a defined variable.
1378 // Note that DefinedVariable can be overloaded
1379
1380 ctemp = chaine;
1381 ctemp.ReplaceAll(escapedSlash, slash);
1382 Int_t action;
1383 k = DefinedVariable(ctemp,action);
1384 if (k==-3) {
1385 // Error message already issued
1386 err = 1;
1387 } else if (k==-2) {
1388 err = 31;
1389 chaine_error = ctemp;
1390 } else if ( k >= 0 ) {
1391 fExpr[fNoper] = ctemp;
1392 actionCode = action;
1393 actionParam = k;
1394 SetAction(fNoper,actionCode,actionParam);
1395 if (action==kDefinedString) fNstring++;
1396 else if (k <kMAXFOUND && !fAlreadyFound.TestBitNumber(k)) {
1397 fAlreadyFound.SetBitNumber(k);
1398 fNval++;
1399 }
1400 fNoper++;
1401 } else if (chaine(0,1) == "!") {
1402 ctemp = chaine(1,lchain-1);
1403 Analyze(ctemp.Data(),err,offset); if (err) return;
1404 fExpr[fNoper] = "!";
1405 actionCode = kNot;
1406 SetAction(fNoper,actionCode,actionParam);
1407 fNoper++;
1408 if (!CheckOperands(fNoper-1,err)) return;
1409 } else if (chaine(0,1)=="\"" && chaine(chaine.Length()-1,1)=="\"") {
1410 // It is a string !!!
1411 fExpr[fNoper] = chaine(1,chaine.Length()-2);
1412 actionCode = kStringConst;
1413 SetAction(fNoper,actionCode,actionParam);
1414 fNoper++;
1415 } else if (chaine(0,4) == "cos(") {
1416 ctemp = chaine(3,lchain-3);
1417 Analyze(ctemp.Data(),err,offset); if (err) return;
1418 fExpr[fNoper] = "cos";
1419 actionCode = kcos;
1420 SetAction(fNoper,actionCode,actionParam);
1421 fNoper++;
1422 if (!CheckOperands(fNoper-1,err)) return;
1423 } else if (chaine(0,4) == "sin(") {
1424 ctemp = chaine(3,lchain-3);
1425 Analyze(ctemp.Data(),err,offset); if (err) return;
1426 fExpr[fNoper] = "sin";
1427 actionCode = ksin;
1428 SetAction(fNoper,actionCode,actionParam);
1429 fNoper++;
1430 if (!CheckOperands(fNoper-1,err)) return;
1431 } else if (chaine(0,4) == "tan(") {
1432 ctemp = chaine(3,lchain-3);
1433 Analyze(ctemp.Data(),err,offset); if (err) return;
1434 fExpr[fNoper] = "tan";
1435 actionCode = ktan;
1436 SetAction(fNoper,actionCode,actionParam);
1437 fNoper++;
1438 if (!CheckOperands(fNoper-1,err)) return;
1439 } else if (chaine(0,5) == "acos(") {
1440 ctemp = chaine(4,lchain-4);
1441 Analyze(ctemp.Data(),err,offset); if (err) return;
1442 fExpr[fNoper] = "acos";
1443 actionCode = kacos;
1444 SetAction(fNoper,actionCode,actionParam);
1445 fNoper++;
1446 if (!CheckOperands(fNoper-1,err)) return;
1447 } else if (chaine(0,5) == "asin(") {
1448 ctemp = chaine(4,lchain-4);
1449 Analyze(ctemp.Data(),err,offset); if (err) return;
1450 fExpr[fNoper] = "asin";
1451 actionCode = kasin;
1452 SetAction(fNoper,actionCode,actionParam);
1453 fNoper++;
1454 if (!CheckOperands(fNoper-1,err)) return;
1455 } else if (chaine(0,5) == "atan(") {
1456 ctemp = chaine(4,lchain-4);
1457 Analyze(ctemp.Data(),err,offset); if (err) return;
1458 fExpr[fNoper] = "atan";
1459 actionCode = katan;
1460 SetAction(fNoper,actionCode,actionParam);
1461 fNoper++;
1462 if (!CheckOperands(fNoper-1,err)) return;
1463 } else if (chaine(0,5) == "cosh(") {
1464 ctemp = chaine(4,lchain-4);
1465 Analyze(ctemp.Data(),err,offset); if (err) return;
1466 fExpr[fNoper] = "cosh";
1467 actionCode = kcosh;
1468 SetAction(fNoper,actionCode,actionParam);
1469 fNoper++;
1470 if (!CheckOperands(fNoper-1,err)) return;
1471 } else if (chaine(0,5) == "sinh(") {
1472 ctemp = chaine(4,lchain-4);
1473 Analyze(ctemp.Data(),err,offset); if (err) return;
1474 fExpr[fNoper] = "sinh";
1475 actionCode = ksinh;
1476 SetAction(fNoper,actionCode,actionParam);
1477 fNoper++;;
1478 if (!CheckOperands(fNoper-1,err)) return;
1479 } else if (chaine(0,5) == "tanh(") {
1480 ctemp = chaine(4,lchain-4);
1481 Analyze(ctemp.Data(),err,offset); if (err) return;
1482 fExpr[fNoper] = "tanh";
1483 actionCode = ktanh;
1484 SetAction(fNoper,actionCode,actionParam);
1485 fNoper++;;
1486 if (!CheckOperands(fNoper-1,err)) return;
1487 } else if (chaine(0,6) == "acosh(") {
1488 ctemp = chaine(5,lchain-5);
1489 Analyze(ctemp.Data(),err,offset); if (err) return;
1490 fExpr[fNoper] = "acosh";
1491 actionCode = kacosh;
1492 SetAction(fNoper,actionCode,actionParam);
1493 fNoper++;;
1494 if (!CheckOperands(fNoper-1,err)) return;
1495 } else if (chaine(0,6) == "asinh(") {
1496 ctemp = chaine(5,lchain-5);
1497 Analyze(ctemp.Data(),err,offset); if (err) return;
1498 fExpr[fNoper] = "asinh";
1499 actionCode = kasinh;
1500 SetAction(fNoper,actionCode,actionParam);
1501 fNoper++;;
1502 if (!CheckOperands(fNoper-1,err)) return;
1503 } else if (chaine(0,6) == "atanh(") {
1504 ctemp = chaine(5,lchain-5);
1505 Analyze(ctemp.Data(),err,offset); if (err) return;
1506 fExpr[fNoper] = "atanh";
1507 actionCode = katanh;
1508 SetAction(fNoper,actionCode,actionParam);
1509 fNoper++;;
1510 if (!CheckOperands(fNoper-1,err)) return;
1511 } else if (chaine(0,3) == "sq(") {
1512 ctemp = chaine(2,lchain-2);
1513 Analyze(ctemp.Data(),err,offset); if (err) return;
1514 fExpr[fNoper] = "sq";
1515 actionCode = ksq;
1516 SetAction(fNoper,actionCode,actionParam);
1517 fNoper++;;
1518 if (!CheckOperands(fNoper-1,err)) return;
1519 } else if (chaine(0,4) == "log(") {
1520 ctemp = chaine(3,lchain-3);
1521 Analyze(ctemp.Data(),err,offset); if (err) return;
1522 fExpr[fNoper] = "log";
1523 actionCode = klog;
1524 SetAction(fNoper,actionCode,actionParam);
1525 fNoper++;;
1526 if (!CheckOperands(fNoper-1,err)) return;
1527 } else if (chaine(0,6) == "log10(") {
1528 ctemp = chaine(5,lchain-5);
1529 Analyze(ctemp.Data(),err,offset); if (err) return;
1530 fExpr[fNoper] = "log10";
1531 actionCode = klog10;
1532 SetAction(fNoper,actionCode,actionParam);
1533 fNoper++;;
1534 if (!CheckOperands(fNoper-1,err)) return;
1535 } else if (chaine(0,4) == "exp(") {
1536 ctemp = chaine(3,lchain-3);
1537 Analyze(ctemp.Data(),err,offset); if (err) return;
1538 fExpr[fNoper] = "exp";
1539 actionCode = kexp;
1540 SetAction(fNoper,actionCode,actionParam);
1541 fNoper++;;
1542 if (!CheckOperands(fNoper-1,err)) return;
1543 } else if (chaine(0,4) == "abs(") {
1544 ctemp = chaine(3,lchain-3);
1545 Analyze(ctemp.Data(),err,offset); if (err) return;
1546 fExpr[fNoper] = "abs";
1547 actionCode = kabs;
1548 SetAction(fNoper,actionCode,actionParam);
1549 fNoper++;;
1550 if (!CheckOperands(fNoper-1,err)) return;
1551 } else if (chaine(0,5) == "sign(") {
1552 ctemp = chaine(4,lchain-4);
1553 Analyze(ctemp.Data(),err,offset); if (err) return;
1554 fExpr[fNoper] = "sign";
1555 actionCode = ksign;
1556 SetAction(fNoper,actionCode,actionParam);
1557 fNoper++;;
1558 if (!CheckOperands(fNoper-1,err)) return;
1559 } else if (chaine(0,4) == "int(") {
1560 ctemp = chaine(3,lchain-3);
1561 Analyze(ctemp.Data(),err,offset); if (err) return;
1562 fExpr[fNoper] = "int";
1563 actionCode = kint;
1564 SetAction(fNoper,actionCode,actionParam);
1565 fNoper++;;
1566 if (!CheckOperands(fNoper-1,err)) return;
1567 } else if (chaine == "rndm" || chaine(0,5) == "rndm(") {
1568 fExpr[fNoper] = "rndm";
1569 actionCode = krndm;
1570 SetAction(fNoper,actionCode,actionParam);
1571 fNoper++;;
1572 } else if (chaine(0,5) == "sqrt(") {
1573 ctemp = chaine(4,lchain-4);
1574 Analyze(ctemp.Data(),err,offset); if (err) return;
1575 fExpr[fNoper] = "sqrt";
1576 actionCode = ksqrt;
1577 SetAction(fNoper,actionCode,actionParam);
1578 fNoper++;;
1579 if (!CheckOperands(fNoper-1,err)) return;
1580
1581 // Look for an exponential
1582
1583 } else if ( chaine == "expo" || chaine(0,5)=="expo("
1584 || (lchain==5 && chaine(1,4)=="expo")
1585 || (lchain==6 && chaine(2,4)=="expo")
1586 || chaine(1,5)=="expo(" || chaine(2,5)=="expo(" ) {
1587 chaine1ST=chaine;
1588 if (chaine(1,4) == "expo") {
1589 ctemp=chaine(0,1);
1590 if (ctemp=="x") {
1591 inter2=0;
1592 if (fNdim < 1) fNdim = 1; }
1593 else if (ctemp=="y") {
1594 inter2=1;
1595 if (fNdim < 2) fNdim = 2; }
1596 else if (ctemp=="z") {
1597 inter2=2;
1598 if (fNdim < 3) fNdim = 3; }
1599 else if (ctemp=="t") {
1600 inter2=3;
1601 if (fNdim < 4) fNdim = 4; }
1602 else {
1603 err=26; // unknown name;
1604 chaine_error=chaine1ST;
1605 }
1606 chaine=chaine(1,lchain-1);
1607 lchain=chaine.Length();
1608 } else inter2=0;
1609 if (chaine(2,4) == "expo") {
1610 if (chaine(0,2) != "xy") {
1611 err=26; // unknown name
1612 chaine_error=chaine1ST;
1613 }
1614 else {
1615 inter2=5;
1616 if (fNdim < 2) fNdim = 2;
1617 chaine=chaine(2,lchain-2);
1618 lchain=chaine.Length();
1619 }
1620 }
1621 if (lchain == 4) {
1622 if (fNpar>=gMAXPAR) err=7; // too many parameters
1623 if (!err) {
1624 fExpr[fNoper] = chaine1ST;
1625 actionCode = kexpo + inter2;
1626 actionParam = offset;
1627 SetAction(fNoper,actionCode,actionParam);
1628 if (inter2 == 5+offset && fNpar < 3+offset) fNpar = 3+offset;
1629 if (fNpar < 2+offset) fNpar = 2+offset;
1630 if (fNpar>=gMAXPAR) err=7; // too many parameters
1631 if (!err) {
1632 fNoper++;
1633 if (fNdim < 1) fNdim = 1;
1634 if (fNpar == 2) SetNumber(200);
1635 }
1636 }
1637 } else if (chaine(4,1) == "(") {
1638 ctemp = chaine(5,lchain-6);
1639 fExpr[fNoper] = chaine1ST;
1640 for (j=0; j<ctemp.Length(); j++) {
1641 t=ctemp[j];
1642 if (strchr("0123456789",t)==0 && (ctemp(j,1)!="+" || j!=0)) {
1643 err=20;
1644 chaine_error=chaine1ST;
1645 }
1646 }
1647 if (err==0) {
1648 sscanf(ctemp.Data(),"%d",&inter);
1649 if (inter>=0) {
1650 inter += offset;
1651 actionCode = kexpo + inter2;
1652 actionParam = inter;
1653 SetAction(fNoper,actionCode,actionParam);
1654 if (inter2 == 5) inter++;
1655 if (inter+2>fNpar) fNpar = inter+2;
1656 if (fNpar>=gMAXPAR) err=7; // too many parameters
1657 if (!err) fNoper++;
1658 if (fNpar == 2) SetNumber(200);
1659 } else err=20;
1660 } else err = 20; // non integer value for parameter number
1661 } else {
1662 err=26; // unknown name
1663 chaine_error=chaine;
1664 }
1665
1666 // Look for gaus, xgaus,ygaus,xygaus
1667
1668 } else if (chaine=="gaus"
1669 || (lchain==5 && chaine(1,4)=="gaus")
1670 || (lchain==6 && chaine(2,4)=="gaus")
1671 || chaine(0,5)=="gaus(" || chaine(1,5)=="gaus(" || chaine(2,5)=="gaus(") {
1672 chaine1ST=chaine;
1673 if (chaine(1,4) == "gaus") {
1674 ctemp=chaine(0,1);
1675 if (ctemp=="x") {
1676 inter2=0;
1677 if (fNdim < 1) fNdim = 1; }
1678 else if (ctemp=="y") {
1679 inter2=1;
1680 if (fNdim < 2) fNdim = 2; }
1681 else if (ctemp=="z") {
1682 inter2=2;
1683 if (fNdim < 3) fNdim = 3; }
1684 else if (ctemp=="t") {
1685 inter2=3;
1686 if (fNdim < 4) fNdim = 4; }
1687 else {
1688 err=26; // unknown name
1689 chaine_error=chaine1ST;
1690 }
1691 chaine=chaine(1,lchain-1);
1692 lchain=chaine.Length();
1693 } else inter2=0;
1694 if (chaine(2,4) == "gaus") {
1695 if (chaine(0,2) != "xy") {
1696 err=26; // unknown name
1697 chaine_error=chaine1ST;
1698 }
1699 else {
1700 inter2=5;
1701 if (fNdim < 2) fNdim = 2;
1702 chaine=chaine(2,lchain-2);
1703 lchain=chaine.Length();
1704 SetNumber(110); // xygaus
1705 }
1706 }
1707 if (lchain == 4 && err==0) {
1708 if (fNpar>=gMAXPAR) err=7; // too many parameters
1709 if (!err) {
1710 fExpr[fNoper] = chaine1ST;
1711 actionCode = kgaus + inter2;
1712 actionParam = offset;
1713 SetAction(fNoper,actionCode,actionParam);
1714 if (inter2 == 5+offset && fNpar < 5+offset) fNpar = 5+offset;
1715 if (3+offset>fNpar) fNpar = 3+offset;
1716 if (fNpar>=gMAXPAR) err=7; // too many parameters
1717 if (!err) {
1718 fNoper++;
1719 if (fNdim < 1) fNdim = 1;
1720 if (fNpar == 3) SetNumber(100);
1721 }
1722 }
1723 } else if (chaine(4,1) == "(" && err==0) {
1724 ctemp = chaine(5,lchain-6);
1725 fExpr[fNoper] = chaine1ST;
1726 for (j=0; j<ctemp.Length(); j++) {
1727 t=ctemp[j];
1728 if (strchr("0123456789",t)==0 && (ctemp(j,1)!="+" || j!=0)) {
1729 err=20;
1730 chaine_error=chaine1ST;
1731 }
1732 }
1733 if (err==0) {
1734 sscanf(ctemp.Data(),"%d",&inter);
1735 if (inter >= 0) {
1736 inter += offset;
1737 actionCode = kgaus + inter2;
1738 actionParam = inter;
1739 SetAction(fNoper,actionCode,actionParam);
1740 if (inter2 == 5) inter += 2;
1741 if (inter+3>fNpar) fNpar = inter+3;
1742 if (fNpar>=gMAXPAR) err=7; // too many parameters
1743 if (!err) fNoper++;
1744 if(fNpar == 3) SetNumber(100);
1745 } else err = 20; // non integer value for parameter number
1746 }
1747 } else if (err==0) {
1748 err=26; // unknown name
1749 chaine_error=chaine1ST;
1750 }
1751
1752 // Look for landau, xlandau,ylandau,xylandau
1753
1754 } else if (chaine=="landau" || (lchain==7 && chaine(1,6)=="landau")
1755 || (lchain==8 && chaine(2,6)=="landau")
1756 || chaine(0,7)=="landau(" || chaine(1,7)=="landau(" || chaine(2,7)=="landau(") {
1757 chaine1ST=chaine;
1758 if (chaine(1,6) == "landau") {
1759 ctemp=chaine(0,1);
1760 if (ctemp=="x") {
1761 inter2=0;
1762 if (fNdim < 1) fNdim = 1; }
1763 else if (ctemp=="y") {
1764 inter2=1;
1765 if (fNdim < 2) fNdim = 2; }
1766 else if (ctemp=="z") {
1767 inter2=2;
1768 if (fNdim < 3) fNdim = 3; }
1769 else if (ctemp=="t") {
1770 inter2=3;
1771 if (fNdim < 4) fNdim = 4; }
1772 else {
1773 err=26; // unknown name
1774 chaine_error=chaine1ST;
1775 }
1776 chaine=chaine(1,lchain-1);
1777 lchain=chaine.Length();
1778 } else inter2=0;
1779 if (chaine(2,6) == "landau") {
1780 if (chaine(0,2) != "xy") {
1781 err=26; // unknown name
1782 chaine_error=chaine1ST;
1783 }
1784 else {
1785 inter2=5;
1786 if (fNdim < 2) fNdim = 2;
1787 chaine=chaine(2,lchain-2);
1788 lchain=chaine.Length();
1789 SetNumber(410);
1790 }
1791 }
1792 if (lchain == 6 && err==0) {
1793 if (fNpar>=gMAXPAR) err=7; // too many parameters
1794 if (!err) {
1795 fExpr[fNoper] = chaine1ST;
1796 actionCode = klandau + inter2;
1797 actionParam = offset;
1798 SetAction(fNoper,actionCode,actionParam);
1799 if (inter2 == 5+offset && fNpar < 5+offset) fNpar = 5+offset;
1800 if (3+offset>fNpar) fNpar = 3+offset;
1801 if (fNpar>=gMAXPAR) err=7; // too many parameters
1802 if (!err) {
1803 fNoper++;
1804 if (fNdim < 1) fNdim = 1;
1805 if (fNpar == 3) SetNumber(400);
1806 }
1807 }
1808 } else if (chaine(6,1) == "(" && err==0) {
1809 ctemp = chaine(7,lchain-8);
1810 fExpr[fNoper] = chaine1ST;
1811 for (j=0; j<ctemp.Length(); j++) {
1812 t=ctemp[j];
1813 if (strchr("0123456789",t)==0 && (ctemp(j,1)!="+" || j!=0)) {
1814 err=20;
1815 chaine_error=chaine1ST;
1816 }
1817 }
1818 if (err==0) {
1819 sscanf(ctemp.Data(),"%d",&inter);
1820 if (inter >= 0) {
1821 inter += offset;
1822 actionCode = klandau + inter2;
1823 actionParam = inter;
1824 SetAction(fNoper,actionCode,actionParam);
1825 if (inter2 == 5) inter += 2;
1826 if (inter+3>fNpar) fNpar = inter+3;
1827 if (fNpar>=gMAXPAR) err=7; // too many parameters
1828 if (!err) fNoper++;
1829 if (fNpar == 3) SetNumber(400);
1830 } else err = 20; // non integer value for parameter number
1831 }
1832 } else if (err==0) {
1833 err=26; // unknown name
1834 chaine_error=chaine1ST;
1835 }
1836
1837 // Look for a polynomial
1838
1839 } else if (chaine(0,3) == "pol" || chaine(1,3) == "pol") {
1840 chaine1ST=chaine;
1841 if (chaine(1,3) == "pol") {
1842 ctemp=chaine(0,1);
1843 if (ctemp=="x") {
1844 inter2=1;
1845 if (fNdim < 1) fNdim = 1; }
1846 else if (ctemp=="y") {
1847 inter2=2;
1848 if (fNdim < 2) fNdim = 2; }
1849 else if (ctemp=="z") {
1850 inter2=3;
1851 if (fNdim < 3) fNdim = 3; }
1852 else if (ctemp=="t") {
1853 inter2=4;
1854 if (fNdim < 4) fNdim = 4; }
1855 else {
1856 err=26; // unknown name;
1857 chaine_error=chaine1ST;
1858 }
1859 chaine=chaine(1,lchain-1);
1860 lchain=chaine.Length();
1861 } else inter2=1;
1862 if (chaine(lchain-1,1) == ")") {
1863 nomb = 0;
1864 for (j=3;j<lchain;j++) if (chaine(j,1)=="(" && nomb == 0) nomb = j;
1865 if (nomb == 3) err = 23; // degree of polynomial not specified
1866 if (nomb == 0) err = 40; // '(' is expected
1867 ctemp = chaine(nomb+1,lchain-nomb-2);
1868 for (j=0; j<ctemp.Length(); j++) {
1869 t=ctemp[j];
1870 if (strchr("0123456789",t)==0 && (ctemp(j,1)!="+" || j!=0)) {
1871 err=20;
1872 chaine_error=chaine1ST;
1873 }
1874 }
1875 if (!err) {
1876 sscanf(ctemp.Data(),"%d",&inter);
1877 if (inter < 0) err = 20;
1878 }
1879 }
1880 else {
1881 nomb = lchain;
1882 inter = 0;
1883 }
1884 if (!err) {
1885 inter--;
1886 ctemp = chaine(3,nomb-3);
1887 if (sscanf(ctemp.Data(),"%d",&n) > 0) {
1888 if (n < 0 ) err = 24; //Degree of polynomial must be positive
1889 if (n >= 20) err = 25; //Degree of polynomial must be less than 20
1890 } else err = 20;
1891 }
1892 if (!err) {
1893 fExpr[fNoper] = chaine1ST;
1894 actionCode = kpol+(inter2-1);
1895 actionParam = n*100+inter+2;
1896 SetAction(fNoper,actionCode,actionParam);
1897 if (inter+n+1>=fNpar) fNpar = inter + n + 2;
1898 if (fNpar>=gMAXPAR) err=7; // too many parameters
1899 if (!err) {
1900 fNoper++;
1901 if (fNdim < 1) fNdim = 1;
1902 SetNumber(300+n);
1903 }
1904 }
1905
1906 // Look for pow,atan2,etc
1907
1908 } else if (chaine(0,4) == "pow(") {
1909 compt = 4; nomb = 0; virgule = 0; nest=0;
1910 while(compt != lchain) {
1911 compt++;
1912 if (chaine(compt-1,1) == "(") nest++;
1913 else if (chaine(compt-1,1) == ")") nest--;
1914 else if (chaine(compt-1,1) == "," && nest==0) {
1915 nomb++;
1916 if (nomb == 1 && virgule == 0) virgule = compt;
1917 }
1918 }
1919 if (nomb != 1) err = 22; // There are plus or minus than 2 arguments for pow
1920 else {
1921 ctemp = chaine(4,virgule-5);
1922 Analyze(ctemp.Data(),err,offset); if (err) return;
1923 UInt_t leftopr = fNoper-1;
1924 ctemp = chaine(virgule,lchain-virgule-1);
1925 Analyze(ctemp.Data(),err,offset); if (err) return;
1926 fExpr[fNoper] = "^";
1927 actionCode = kpow;
1928 SetAction(fNoper,actionCode,actionParam);
1929 fNoper++;
1930 if (!CheckOperands(leftopr,fNoper-1,err)) return;
1931 }
1932 } else if (chaine(0,7) == "strstr(") {
1933 compt = 7; nomb = 0; virgule = 0; nest=0;
1934 inString = false;
1935 while(compt != lchain) {
1936 compt++;
1937 if (chaine(compt-1,1) == "\"") {
1938 inString = !inString;
1939 } else if (!inString) {
1940 if (chaine(compt-1,1) == "(") nest++;
1941 else if (chaine(compt-1,1) == ")") nest--;
1942 else if (chaine(compt-1,1) == "," && nest==0) {
1943 nomb++;
1944 if (nomb == 1 && virgule == 0) virgule = compt;
1945 }
1946 }
1947 }
1948 if (nomb != 1) err = 28; // There are plus or minus than 2 arguments for strstr
1949 else {
1950 ctemp = chaine(7,virgule-8);
1951 Analyze(ctemp.Data(),err,offset); if (err) return;
1952 Int_t optloc = fNoper-1;
1953
1954 ctemp = chaine(virgule,lchain-virgule-1);
1955 Analyze(ctemp.Data(),err,offset); if (err) return;
1956 fExpr[fNoper] = "strstr";
1957 actionCode = kstrstr;
1958 SetAction(fNoper,actionCode,actionParam);
1959 fNoper++;
1960
1961 if ( !IsString(optloc) || !IsString(fNoper-2) ) {
1962 err = 46;
1963 chaine_error = "strstr";
1964 }
1965 }
1966 } else if (chaine(0,4) == "min(") {
1967 compt = 4; nomb = 0; virgule = 0; nest=0;
1968 while(compt != lchain) {
1969 compt++;
1970 if (chaine(compt-1,1) == "(") nest++;
1971 else if (chaine(compt-1,1) == ")") nest--;
1972 else if (chaine(compt-1,1) == "," && nest==0) {
1973 nomb++;
1974 if (nomb == 1 && virgule == 0) virgule = compt;
1975 }
1976 }
1977 if (nomb != 1) {
1978 err = 44; // There are plus or minus than 2 arguments for min
1979 err_hint = 3;
1980 }
1981 else {
1982 ctemp = chaine(4,virgule-5);
1983 Analyze(ctemp.Data(),err,offset); if (err) return;
1984 UInt_t leftopr = fNoper-1;
1985 ctemp = chaine(virgule,lchain-virgule-1);
1986 Analyze(ctemp.Data(),err,offset); if (err) return;
1987 fExpr[fNoper] = "min";
1988 actionCode = kmin;
1989 SetAction(fNoper,actionCode,actionParam);
1990 fNoper++;
1991 if (!CheckOperands(leftopr,fNoper-1,err)) return;
1992 }
1993 } else if (chaine(0,4) == "max(") {
1994 compt = 4; nomb = 0; virgule = 0; nest=0;
1995 while(compt != lchain) {
1996 compt++;
1997 if (chaine(compt-1,1) == "(") nest++;
1998 else if (chaine(compt-1,1) == ")") nest--;
1999 else if (chaine(compt-1,1) == "," && nest==0) {
2000 nomb++;
2001 if (nomb == 1 && virgule == 0) virgule = compt;
2002 }
2003 }
2004 if (nomb != 1) {
2005 err = 44; // There are plus or minus than 2 arguments for min
2006 err_hint = 3;
2007 }
2008 else {
2009 ctemp = chaine(4,virgule-5);
2010 Analyze(ctemp.Data(),err,offset); if (err) return;
2011 UInt_t leftopr = fNoper-1;
2012 ctemp = chaine(virgule,lchain-virgule-1);
2013 Analyze(ctemp.Data(),err,offset); if (err) return;
2014 fExpr[fNoper] = "max";
2015 actionCode = kmax;
2016 SetAction(fNoper,actionCode,actionParam);
2017 fNoper++;
2018 if (!CheckOperands(leftopr,fNoper-1,err)) return;
2019 }
2020
2021 } else if (chaine(0,6) == "atan2(") {
2022 compt = 6; nomb = 0; virgule = 0; nest=0;
2023 while(compt != lchain) {
2024 compt++;
2025 if (chaine(compt-1,1) == "(") nest++;
2026 else if (chaine(compt-1,1) == ")") nest--;
2027 else if (chaine(compt-1,1) == "," && nest==0) {
2028 nomb++;
2029 if (nomb == 1 && virgule == 0) virgule = compt;
2030 }
2031 }
2032 if (nomb != 1) err = 21; //{ There are plus or minus than 2 arguments for atan2
2033 else {
2034 ctemp = chaine(6,virgule-7);
2035 Analyze(ctemp.Data(),err,offset); if (err) return;
2036 UInt_t leftopr = fNoper-1;
2037 ctemp = chaine(virgule,lchain-virgule-1);
2038 Analyze(ctemp.Data(),err,offset); if (err) return;
2039 fExpr[fNoper] = "atan2";
2040 actionCode = katan2;
2041 SetAction(fNoper,actionCode,actionParam);
2042 fNoper++;
2043 if (!CheckOperands(leftopr,fNoper-1,err)) return;
2044 }
2045 } else if (chaine(0,5) == "fmod(") {
2046 compt = 5; nomb = 0; virgule = 0; nest=0;
2047 while(compt != lchain) {
2048 compt++;
2049 if (chaine(compt-1,1) == "(") nest++;
2050 else if (chaine(compt-1,1) == ")") nest--;
2051 else if (chaine(compt-1,1) == "," && nest==0) {
2052 nomb++;
2053 if (nomb == 1 && virgule == 0) virgule = compt;
2054 }
2055 }
2056 if (nomb != 1) {
2057 err = 44; // There are plus or minus than 2 arguments for fmod
2058 err_hint = 4;
2059 }
2060 else {
2061 ctemp = chaine(5,virgule-6);
2062 Analyze(ctemp.Data(),err,offset); if (err) return;
2063 UInt_t leftopr = fNoper-1;
2064 ctemp = chaine(virgule,lchain-virgule-1);
2065 Analyze(ctemp.Data(),err,offset); if (err) return;
2066 fExpr[fNoper] = "fmod";
2067 actionCode = kfmod;
2068 SetAction(fNoper,actionCode,actionParam);
2069 fNoper++;
2070 if (!CheckOperands(leftopr,fNoper-1,err)) return;
2071 }
2072 } else if (AnalyzeFunction(chaine,err,offset) || err) { // The '||err' is to grab an error coming from AnalyzeFunction
2073 if (err) {
2074 chaine_error = chaine;
2075 } else {
2076 // We have a function call. Note that all the work was already,
2077 // eventually done in AnalyzeFunction
2078 //fprintf(stderr,"We found a foreign function in %s\n",chaine.Data());
2079 }
2080 } else if (chaine(0,1) == "[" && chaine(lchain-1,1) == "]") {
2081 fExpr[fNoper] = chaine;
2082 fNoper++;
2083 ctemp = chaine(1,lchain-2);
2084 for (j=0; j<ctemp.Length(); j++) {
2085 t=ctemp[j];
2086 if (strchr("0123456789",t)==0 && (ctemp(j,1)!="+" || j!=0)) {
2087 err=20;
2088 chaine_error=chaine1ST; // le numero ? de par[?] n'est pas un entier }
2089 }
2090 }
2091 if (!err) {
2092 sscanf(ctemp.Data(),"%d",&valeur);
2093 actionCode = kParameter;
2094 actionParam = offset + valeur;
2095 SetAction(fNoper-1, actionCode, actionParam);
2096 fExpr[fNoper-1] = "[";
2097 fExpr[fNoper-1] = (fExpr[fNoper-1] + (long int)(valeur+offset)) + "]";
2098 }
2099 } else if (chaine == "pi") {
2100 fExpr[fNoper] = "pi";
2101 actionCode = kpi;
2102 SetAction(fNoper,actionCode,actionParam);
2103 fNoper++;
2104 }
2105 else {
2106
2107 // None of the above.
2108
2109 err = 30;
2110 }
2111 }
2112 }
2113 }
2114 }
2115 }
2116 }
2117
2118 // Overflows
2119 if (fNoper>=gMAXOP) err=6; // too many operators
2120
2121 }
2122
2123 // errors!
2124 if (err>1) {
2125 TString er = "";
2126 chaine_error = "\""+chaine_error+"\"";
2127 switch(err) {
2128 case 2 : er = " Invalid Floating Point Operation"; break;
2129 case 4 : er = " Empty String"; break;
2130 case 5 : er = " Invalid Syntax " + chaine_error; break;
2131 case 6 : er = " Too many operators !"; break;
2132 case 7 : er = " Too many parameters !"; break;
2133 case 10 : er = " z specified but not x and y"; break;
2134 case 11 : er = " z and y specified but not x"; break;
2135 case 12 : er = " y specified but not x"; break;
2136 case 13 : er = " z and x specified but not y"; break;
2137 case 20 : er = " Non integer value for parameter number : " + chaine_error; break;
2138 case 21 : er = " ATAN2 requires two arguments"; break;
2139 case 22 : er = " POW requires two arguments"; break;
2140 case 23 : er = " Degree of polynomial not specified"; break;
2141 case 24 : er = " Degree of polynomial must be positive"; break;
2142 case 25 : er = " Degree of polynomial must be less than 20"; break;
2143 case 26 : er = " Unknown name : " + chaine_error; break;
2144 case 27 : er = " Too many constants in expression"; break;
2145 case 28 : er = " strstr requires two arguments"; break;
2146 case 29 : er = " TFormula can only call interpreted and compiled functions that return a numerical type: " + chaine_error; break;
2147 case 30 : er = " Bad numerical expression : " + chaine_error; break;
2148 case 31 : er = " Part of the Variable " + chaine_error; er += " exists but some of it is not accessible or useable"; break;
2149 case 40 : er = " '(' is expected"; break;
2150 case 41 : er = " ')' is expected"; break;
2151 case 42 : er = " '[' is expected"; break;
2152 case 43 : er = " ']' is expected"; break;
2153 case 44 : er = " The function '" + chaine(0,err_hint) + "' requires two arguments."; break;
2154 case 45 : er = "The operator " + chaine_error + " requires a numerical operand."; break;
2155 case 46 : er = "Both operands of the operator " + chaine_error + " have to be either numbers or strings."; break;
2156 case 47 : er = chaine_error + " requires 2 string arguments"; break;
2157 }
2158 Error("Compile", "%s", er.Data());
2159 err=1;
2160 }
2161
2162}
2163
2164////////////////////////////////////////////////////////////////////////////////
2165/// Check whether the operand at 'oper-1' is compatible with the operation
2166/// at 'oper'.
2167
2168Bool_t TFormula::CheckOperands(Int_t oper, Int_t &err)
2169{
2170 if ( IsString(oper-1) && !StringToNumber(oper-1) ) {
2171 Error("Compile","\"%s\" requires a numerical operand.",fExpr[oper].Data());
2172 err = 45;
2173 return kFALSE;
2174 }
2175 return kTRUE;
2176}
2177
2178////////////////////////////////////////////////////////////////////////////////
2179/// Check whether the operands at 'leftoper' and 'oper-1' are compatible with
2180/// the operation at 'oper'.
2181
2182Bool_t TFormula::CheckOperands(Int_t leftoper, Int_t oper, Int_t &err)
2183{
2184 if ( IsString(oper-1) || IsString(leftoper) ) {
2185 if (IsString(oper-1) && StringToNumber(oper-1)) {
2186 return kTRUE;
2187 }
2188 if (IsString(leftoper) && StringToNumber(leftoper)) {
2189 return kTRUE;
2190 }
2191 Error("Compile","\"%s\" requires two numerical operands.",fExpr[oper].Data());
2192 err = 46;
2193 return kFALSE;
2194 }
2195 return kTRUE;
2196}
2197
2198////////////////////////////////////////////////////////////////////////////////
2199/// Try to 'demote' a string into an array bytes. If this is not possible,
2200/// return false.
2201
2202Bool_t TFormula::StringToNumber(Int_t /* code */)
2203{
2204 // In TFormula proper, we can not handle array of bytes ...
2205 return kFALSE;
2206}
2207
2208////////////////////////////////////////////////////////////////////////////////
2209/// Resets the objects.
2210///
2211/// Resets the object to its state before compilation.
2212
2213void TFormula::Clear(Option_t * /*option*/ )
2214{
2215 ClearFormula();
2216}
2217
2218////////////////////////////////////////////////////////////////////////////////
2219/// Resets the objects.
2220///
2221/// Resets the object to its state before compilation.
2222
2223void TFormula::ClearFormula(Option_t * /*option*/ )
2224{
2225 fNdim = 0;
2226 fNpar = 0;
2227 fNoper = 0;
2228 fNconst = 0;
2229 fNumber = 0;
2230 fNstring= 0;
2231 fNval = 0;
2232
2233 if (fExpr) { delete [] fExpr; fExpr = 0;}
2234 if (fNames) { delete [] fNames; fNames = 0;}
2235 if (fOper) { delete [] fOper; fOper = 0;}
2236 if (fConst) { delete [] fConst; fConst = 0;}
2237 if (fParams) { delete [] fParams; fParams = 0;}
2238 fFunctions.Delete();
2239 fLinearParts.Delete();
2240 //
2241 //MI change
2242 if (fPredefined) { delete [] fPredefined; fPredefined = 0;}
2243 if (fOperOffset) { delete [] fOperOffset; fOperOffset = 0;}
2244 if (fExprOptimized) { delete [] fExprOptimized; fExprOptimized = 0;}
2245 if (fOperOptimized) { delete [] fOperOptimized; fOperOptimized = 0;}
2246 // should we also remove the object from the list?
2247 // gROOT->GetListOfFunctions()->Remove(this);
2248 // if we don't, what happens if it fails the new compilation?
2249}
2250
2251namespace {
2252 template <class T>
2253 inline static void ResizeArrayIfAllocated(T*& oldArray, int newSize){
2254
2255 // Don't do anything in this case.
2256 if (!oldArray || newSize <=0) return;
2257
2258 T* newArray = new T[newSize];
2259 std::copy(oldArray, oldArray+newSize, newArray);
2260 delete [] oldArray;
2261 oldArray = newArray;
2262 }
2263}
2264
2265////////////////////////////////////////////////////////////////////////////////
2266/// Compile expression already stored in fTitle.
2267///
2268/// Loop on all subexpressions of formula stored in fTitle
2269///
2270/// If you overload this member function, you also HAVE TO
2271/// never call the constructor:
2272///
2273/// ~~~ {.cpp}
2274/// TFormula::TFormula(const char *name,const char *expression)
2275/// ~~~
2276///
2277/// and write your own constructor
2278///
2279/// ~~~ {.cpp}
2280/// MyClass::MyClass(const char *name,const char *expression) : TFormula()
2281/// ~~~
2282///
2283/// which has to call the TFormula default constructor and whose implementation
2284/// should be similar to the implementation of the normal TFormula constructor
2285///
2286/// This is necessary because the normal TFormula constructor call indirectly
2287/// the virtual member functions Analyze, DefaultString, DefaultValue
2288/// and DefaultVariable.
2289///
2290/// \image html TFormula_compile.png
2291
2292Int_t TFormula::Compile(const char *expression)
2293{
2294 Int_t i,j,lc,valeur,err;
2295 TString ctemp;
2296
2297 ClearFormula();
2298
2299 // If expression is not empty, take it, otherwise take the title
2300 if (strlen(expression)) SetTitle(expression);
2301
2302 TString chaine = GetTitle();
2303
2304 if (chaine.Contains(";")) {
2305 char *sctemp = new char[chaine.Length()+1];
2306 strlcpy(sctemp,chaine.Data(),chaine.Length()+1);
2307 char *semicol = (char*)strstr(sctemp,";");
2308 if (semicol) *semicol = 0;
2309 chaine = sctemp;
2310 delete [] sctemp;
2311 }
2312
2313 // if the function is linear, process it and fill the array of linear parts
2314 if (TestBit(kLinear)){
2315 ProcessLinear(chaine);
2316 }
2317
2318 // see static function SetMaxima to change the gMAX.. default values
2319 fExpr = new TString[gMAXOP];
2320 fConst = new Double_t[gMAXCONST];
2321 fParams = new Double_t[gMAXPAR];
2322 fNames = new TString[gMAXPAR];
2323 fOper = new Int_t[gMAXOP];
2324 for (i=0; i<gMAXPAR; i++) {
2325 fParams[i] = 0;
2326 fNames[i] = "";
2327 }
2328 for (i=0; i<gMAXOP; i++) {
2329 fExpr[i] = "";
2330 fOper[i] = 0;
2331 }
2332 for (i=0; i<gMAXCONST; i++)
2333 fConst[i] = 0;
2334
2335 // Substitution of some operators to C++ style
2336 Bool_t inString = false;
2337 for (i=1; i<=chaine.Length(); i++) {
2338 lc =chaine.Length();
2339 if (chaine(i-1,1) == "\"") inString = !inString;
2340 if (inString) continue;
2341 if (chaine(i-1,2) == "**") {
2342 chaine = chaine(0,i-1) + "^" + chaine(i+1,lc-i-1);
2343 i=0;
2344 } else if (chaine(i-1,2) == "++") {
2345 chaine = chaine(0,i) + chaine(i+1,lc-i-1);
2346 i=0;
2347 } else if (chaine(i-1,2) == "+-" || chaine(i-1,2) == "-+") {
2348 chaine = chaine(0,i-1) + "-" + chaine(i+1,lc-i-1);
2349 i=0;
2350 } else if (chaine(i-1,2) == "--") {
2351 chaine = chaine(0,i-1) + "+" + chaine(i+1,lc-i-1);
2352 i=0;
2353 } else if (chaine(i-1,2) == "->") {
2354 chaine = chaine(0,i-1) + "." + chaine(i+1,lc-i-1);
2355 i=0;
2356 } else if (chaine(i-1,1) == "[") {
2357 for (j=1;j<=chaine.Length()-i;j++) {
2358 if (chaine(j+i-1,1) == "]" || j+i > chaine.Length()) break;
2359 }
2360 ctemp = chaine(i,j-1);
2361 valeur=0;
2362 sscanf(ctemp.Data(),"%d",&valeur);
2363 if (valeur >= fNpar) fNpar = valeur+1;
2364 } else if (chaine(i-1,1) == " ") {
2365 chaine = chaine(0,i-1)+chaine(i,lc-i);
2366 i=0;
2367 }
2368 }
2369 err = 0;
2370 Analyze((const char*)chaine,err);
2371
2372 // if no parameters delete arrays fParams and fNames
2373 if (!fNpar) {
2374 delete [] fParams; fParams = 0;
2375 delete [] fNames; fNames = 0;
2376 }
2377
2378 // if no errors, copy local parameters to formula objects
2379 if (!err) {
2380 if (fNdim <= 0) fNdim = 1;
2381 if (chaine.Length() > 4)
2382 {
2383 if ( GetNumber() != 400 &&
2384 GetNumber() != 410 &&
2385 GetNumber() != 110 )
2386 SetNumber(0);
2387 else if ( GetNumber() == 110 && chaine.Length() > 6 )
2388 SetNumber(0);
2389 else if ( GetNumber() == 410 && chaine.Length() > 8 )
2390 SetNumber(0);
2391 }
2392 // if formula is a gaussian, set parameter names
2393 if (GetNumber() == 100) {
2394 SetParName(0,"Constant");
2395 SetParName(1,"Mean");
2396 SetParName(2,"Sigma");
2397 }
2398 // if formula is a 2D gaussian, set parameter names
2399 if (GetNumber() == 110){
2400 SetParName(0,"Constant");
2401 SetParName(1,"MeanX");
2402 SetParName(2,"SigmaX");
2403 SetParName(3,"MeanY");
2404 SetParName(4,"SigmaY");
2405 }
2406 // if formula is an exponential, set parameter names
2407 if (GetNumber() == 200) {
2408 SetParName(0,"Constant");
2409 SetParName(1,"Slope");
2410 }
2411 // if formula is a polynom, set parameter names
2412 if (GetNumber() == 300+fNpar) {
2413 for (i=0;i<fNpar;i++) SetParName(i,Form("p%d",i));
2414 }
2415 // if formula is a landau, set parameter names
2416 if (GetNumber() == 400) {
2417 SetParName(0,"Constant");
2418 SetParName(1,"MPV");
2419 SetParName(2,"Sigma");
2420 }
2421 // if formula is a 2D landau, set parameter names
2422 if (GetNumber() == 410) {
2423 SetParName(0,"Constant");
2424 SetParName(1,"MPVX");
2425 SetParName(2,"SigmaX");
2426 SetParName(3,"MPVY");
2427 SetParName(4,"SigmaY");
2428 }
2429 }
2430
2431 // Here we shrink the arrays allocated like this:
2432 // fExpr = new TString[gMAXOP];
2433 // fConst = new Double_t[gMAXCONST];
2434 // fParams = new Double_t[gMAXPAR];
2435 // fNames = new TString[gMAXPAR];
2436 // fOper = new Int_t[gMAXOP];
2437 // fParams and fNames may be already 0, so we have to check.
2438 if (!err){
2439 ResizeArrayIfAllocated(fExpr, fNoper);
2440 ResizeArrayIfAllocated(fConst, fNconst);
2441 ResizeArrayIfAllocated(fParams, fNpar);
2442 ResizeArrayIfAllocated(fNames, fNpar);
2443 ResizeArrayIfAllocated(fOper, fNoper);
2444 }
2445
2446
2447 if (err) { fNdim = 0; return 1; }
2448 // Convert(5);
2449 //
2450 //MI change
2451 if (!IsA()->GetBaseClass("TTreeFormula")) {
2452 Optimize();
2453 }
2454 //
2455 return 0;
2456}
2457
2458////////////////////////////////////////////////////////////////////////////////
2459/// Copy this formula.
2460
2461void TFormula::Copy(TObject &obj) const
2462{
2463 Int_t i;
2464 ((TFormula&)obj).ClearFormula();
2465 TNamed::Copy(obj);
2466 ((TFormula&)obj).fNdim = fNdim;
2467 ((TFormula&)obj).fNpar = fNpar;
2468 ((TFormula&)obj).fNoper = fNoper;
2469 ((TFormula&)obj).fNconst = fNconst;
2470 ((TFormula&)obj).fNumber = fNumber;
2471 ((TFormula&)obj).fNval = fNval;
2472 ((TFormula&)obj).fExpr = 0;
2473 ((TFormula&)obj).fConst = 0;
2474 ((TFormula&)obj).fParams = 0;
2475 ((TFormula&)obj).fNames = 0;
2476 if (fExpr && fNoper) {
2477 ((TFormula&)obj).fExpr = new TString[fNoper];
2478 for (i=0;i<fNoper;i++) ((TFormula&)obj).fExpr[i] = fExpr[i];
2479 }
2480 if (fOper && fNoper) {
2481 ((TFormula&)obj).fOper = new Int_t[fNoper];
2482 for (i=0;i<fNoper;i++) ((TFormula&)obj).fOper[i] = fOper[i];
2483 }
2484 if (fConst && fNconst) {
2485 ((TFormula&)obj).fConst = new Double_t[fNconst];
2486 for (i=0;i<fNconst;i++) ((TFormula&)obj).fConst[i] = fConst[i];
2487 }
2488 if (fParams && fNpar) {
2489 ((TFormula&)obj).fParams = new Double_t[fNpar];
2490 for (i=0;i<fNpar;i++) ((TFormula&)obj).fParams[i] = fParams[i];
2491 }
2492 if (fNames && fNpar) {
2493 ((TFormula&)obj).fNames = new TString[fNpar];
2494 for (i=0;i<fNpar;i++) ((TFormula&)obj).fNames[i] = fNames[i];
2495 }
2496
2497 TIter next(&fFunctions);
2498 TObject *fobj;
2499 while ( (fobj = next()) ) {
2500 ((TFormula&)obj).fFunctions.Add( fobj->Clone() );
2501 }
2502 //
2503 // MI change
2504 //
2505 //
2506 if (fNoper) {
2507 if(fExprOptimized) {
2508 ((TFormula&)obj).fExprOptimized = new TString[fNoper];
2509 for (i=0;i<fNoper;i++) ((TFormula&)obj).fExprOptimized[i] = fExprOptimized[i];
2510 }
2511 if (fOperOptimized) {
2512 ((TFormula&)obj).fOperOptimized = new Int_t[fNoper];
2513 for (i=0;i<fNoper;i++) ((TFormula&)obj).fOperOptimized[i] = fOperOptimized[i];
2514 }
2515 if (fPredefined) {
2516 ((TFormula&)obj).fPredefined = new ROOT::v5::TFormulaPrimitive*[fNoper];
2517 for (i=0;i<fNoper;i++) {((TFormula&)obj).fPredefined[i] = fPredefined[i];}
2518 }
2519 if (fOperOffset) {
2520 ((TFormula&)obj).fOperOffset = new TOperOffset[fNoper];
2521 for (i=0;i<fNoper;i++) {((TFormula&)obj).fOperOffset[i] = fOperOffset[i];}
2522 }
2523 }
2524 ((TFormula&)obj).fNOperOptimized = fNOperOptimized;
2525 ((TFormula&)obj).fOptimal = fOptimal;
2526
2527}
2528
2529////////////////////////////////////////////////////////////////////////////////
2530/// Return address of string corresponding to special code.
2531///
2532/// This member function is inactive in the TFormula class.
2533/// It may be redefined in derived classes.
2534///
2535/// If you overload this member function, you also HAVE TO
2536/// never call the constructor:
2537///
2538/// ~~~ {.cpp}
2539/// TFormula::TFormula(const char *name,const char *expression)
2540/// ~~~
2541///
2542/// and write your own constructor
2543///
2544/// ~~~ {.cpp}
2545/// MyClass::MyClass(const char *name,const char *expression) : TFormula()
2546/// ~~~
2547///
2548/// which has to call the TFormula default constructor and whose implementation
2549/// should be similar to the implementation of the normal TFormula constructor
2550///
2551/// This is necessary because the normal TFormula constructor call indirectly
2552/// the virtual member functions Analyze, DefaultString, DefaultValue
2553/// and DefaultVariable.
2554
2555char *TFormula::DefinedString(Int_t)
2556{
2557 return 0;
2558}
2559
2560////////////////////////////////////////////////////////////////////////////////
2561/// Return value corresponding to special code.
2562///
2563/// This member function is inactive in the TFormula class.
2564/// It may be redefined in derived classes.
2565///
2566/// If you overload this member function, you also HAVE TO
2567/// never call the constructor:
2568///
2569/// ~~~ {.cpp}
2570/// TFormula::TFormula(const char *name,const char *expression)
2571/// ~~~
2572///
2573/// and write your own constructor
2574///
2575/// ~~~ {.cpp}
2576/// MyClass::MyClass(const char *name,const char *expression) : TFormula()
2577/// ~~~
2578///
2579/// which has to call the TFormula default constructor and whose implementation
2580/// should be similar to the implementation of the normal TFormula constructor
2581///
2582/// This is necessary because the normal TFormula constructor call indirectly
2583/// the virtual member functions Analyze, DefaultString, DefaultValue
2584/// and DefaultVariable.
2585
2586Double_t TFormula::DefinedValue(Int_t)
2587{
2588 return 0;
2589}
2590
2591////////////////////////////////////////////////////////////////////////////////
2592/// Check if expression is in the list of defined variables.
2593///
2594/// This member function can be overloaded in derived classes
2595///
2596/// If you overload this member function, you also HAVE TO
2597/// never call the constructor:
2598///
2599/// ~~~ {.cpp}
2600/// TFormula::TFormula(const char *name,const char *expression)
2601/// ~~~
2602///
2603/// and write your own constructor
2604///
2605/// ~~~ {.cpp}
2606/// MyClass::MyClass(const char *name,const char *expression) : TFormula()
2607/// ~~~
2608///
2609/// which has to call the TFormula default constructor and whose implementation
2610/// should be similar to the implementation of the normal TFormula constructor
2611///
2612/// This is necessary because the normal TFormula constructor call indirectly
2613/// the virtual member functions Analyze, DefaultString, DefaultValue
2614/// and DefaultVariable.
2615///
2616/// The expected returns values are
2617/// - -2 : the name has been recognized but won't be usable
2618/// - -1 : the name has not been recognized
2619/// - >=0 : the name has been recognized, return the action parameter.
2620
2621Int_t TFormula::DefinedVariable(TString &chaine,Int_t &action)
2622{
2623 action = kVariable;
2624 if (chaine == "x") {
2625 if (fNdim < 1) fNdim = 1;
2626 return 0;
2627 } else if (chaine == "y") {
2628 if (fNdim < 2) fNdim = 2;
2629 return 1;
2630 } else if (chaine == "z") {
2631 if (fNdim < 3) fNdim = 3;
2632 return 2;
2633 } else if (chaine == "t") {
2634 if (fNdim < 4) fNdim = 4;
2635 return 3;
2636 }
2637 // MI change
2638 // extended defined variable (MI)
2639 //
2640 if (chaine.Data()[0]=='x'){
2641 if (chaine.Data()[1]=='[' && chaine.Data()[3]==']'){
2642 const char ch0 = '0';
2643 Int_t dim = chaine.Data()[2]-ch0;
2644 if (dim<0) return -1;
2645 if (dim>9) return -1;
2646 if (fNdim<=dim) fNdim = dim+1;
2647 return dim;
2648 }
2649 if (chaine.Data()[1]=='[' && chaine.Data()[4]==']'){
2650 const char ch0 = '0';
2651 Int_t dim = (chaine.Data()[2]-ch0)*10+(chaine.Data()[3]-ch0);
2652 if (dim<0) return -1;
2653 if (dim>99) return -1;
2654 if (fNdim<=dim) fNdim = dim+1;
2655 return dim;
2656 }
2657 }
2658 return -1;
2659}
2660
2661////////////////////////////////////////////////////////////////////////////////
2662/// Evaluate this formula.
2663///
2664/// The current value of variables x,y,z,t is passed through x, y, z and t.
2665/// The parameters used will be the ones in the array params if params is given
2666/// otherwise parameters will be taken from the stored data members fParams
2667
2668Double_t TFormula::Eval(Double_t x, Double_t y, Double_t z, Double_t t) const
2669{
2670 Double_t xx[4];
2671 xx[0] = x;
2672 xx[1] = y;
2673 xx[2] = z;
2674 xx[3] = t;
2675 return ((TFormula*)this)->EvalPar(xx);
2676}
2677
2678////////////////////////////////////////////////////////////////////////////////
2679/// Evaluate this formula.
2680///
2681/// The current value of variables x,y,z,t is passed through the pointer x.
2682/// The parameters used will be the ones in the array params if params is given
2683/// otherwise parameters will be taken from the stored data members fParams
2684///
2685/// \image html TFormula_eval.png
2686
2687Double_t TFormula::EvalParOld(const Double_t *x, const Double_t *uparams)
2688{
2689 Int_t i,j;
2690 // coverity[uninit] the tab value of tab is guaranteed to be set properly by the control flow.
2691 Double_t tab[kMAXFOUND];
2692 const char *stringStack[gMAXSTRINGFOUND] = {0};
2693 Double_t param_calc[kMAXFOUND];
2694 char *string_calc[gMAXSTRINGFOUND] = {0};
2695 Int_t precalculated = 0;
2696 Int_t precalculated_str = 0;
2697 Double_t *params;
2698
2699 if (uparams) {
2700 params = const_cast<Double_t*>(uparams);
2701 } else {
2702 params = fParams;
2703 }
2704 UInt_t pos = 0;
2705 UInt_t strpos = 0;
2706
2707 for (i=0; i<fNoper; ++i) {
2708
2709 const int oper = fOper[i];
2710 const int opcode = oper >> kTFOperShift;
2711
2712 switch(opcode) {
2713
2714 case kParameter : { pos++; tab[pos-1] = params[ oper & kTFOperMask ]; continue; }
2715 case kConstant : { pos++; tab[pos-1] = fConst[ oper & kTFOperMask ]; continue; }
2716 case kVariable : { pos++; tab[pos-1] = x[ oper & kTFOperMask ]; continue; }
2717 case kStringConst: { strpos++; stringStack[strpos-1] = (char*)fExpr[i].Data(); pos++; tab[pos-1] = 0; continue; }
2718
2719 case kAdd : pos--; tab[pos-1] += tab[pos]; continue;
2720 case kSubstract : pos--; tab[pos-1] -= tab[pos]; continue;
2721 case kMultiply : pos--; tab[pos-1] *= tab[pos]; continue;
2722 case kDivide : pos--; if (tab[pos] == 0) tab[pos-1] = 0; // division by 0
2723 else tab[pos-1] /= tab[pos];
2724 continue;
2725 case kModulo : {pos--;
2726 Long64_t int1((Long64_t)tab[pos-1]);
2727 Long64_t int2((Long64_t)tab[pos]);
2728 tab[pos-1] = Double_t(int1%int2);
2729 continue;}
2730
2731 case kcos : tab[pos-1] = TMath::Cos(tab[pos-1]); continue;
2732 case ksin : tab[pos-1] = TMath::Sin(tab[pos-1]); continue;
2733 case ktan : if (TMath::Cos(tab[pos-1]) == 0) {tab[pos-1] = 0;} // { tangente indeterminee }
2734 else tab[pos-1] = TMath::Tan(tab[pos-1]);
2735 continue;
2736 case kacos : if (TMath::Abs(tab[pos-1]) > 1) {tab[pos-1] = 0;} // indetermination
2737 else tab[pos-1] = TMath::ACos(tab[pos-1]);
2738 continue;
2739 case kasin : if (TMath::Abs(tab[pos-1]) > 1) {tab[pos-1] = 0;} // indetermination
2740 else tab[pos-1] = TMath::ASin(tab[pos-1]);
2741 continue;
2742 case katan : tab[pos-1] = TMath::ATan(tab[pos-1]); continue;
2743 case kcosh : tab[pos-1] = TMath::CosH(tab[pos-1]); continue;
2744 case ksinh : tab[pos-1] = TMath::SinH(tab[pos-1]); continue;
2745 case ktanh : if (TMath::CosH(tab[pos-1]) == 0) {tab[pos-1] = 0;} // { tangente indeterminee }
2746 else tab[pos-1] = TMath::TanH(tab[pos-1]);
2747 continue;
2748 case kacosh: if (tab[pos-1] < 1) {tab[pos-1] = 0;} // indetermination
2749 else tab[pos-1] = TMath::ACosH(tab[pos-1]);
2750 continue;
2751 case kasinh: tab[pos-1] = TMath::ASinH(tab[pos-1]); continue;
2752 case katanh: if (TMath::Abs(tab[pos-1]) > 1) {tab[pos-1] = 0;} // indetermination
2753 else tab[pos-1] = TMath::ATanH(tab[pos-1]);
2754 continue;
2755 case katan2: pos--; tab[pos-1] = TMath::ATan2(tab[pos-1],tab[pos]); continue;
2756
2757 case kfmod : pos--; tab[pos-1] = fmod(tab[pos-1],tab[pos]); continue;
2758 case kpow : pos--; tab[pos-1] = TMath::Power(tab[pos-1],tab[pos]); continue;
2759 case ksq : tab[pos-1] = tab[pos-1]*tab[pos-1]; continue;
2760 case ksqrt : tab[pos-1] = TMath::Sqrt(TMath::Abs(tab[pos-1])); continue;
2761
2762 case kstrstr : strpos -= 2; pos-=2; pos++;
2763 if (strstr(stringStack[strpos],stringStack[strpos+1])) tab[pos-1]=1;
2764 else tab[pos-1]=0;
2765 continue;
2766
2767 case kmin : pos--; tab[pos-1] = TMath::Min(tab[pos-1],tab[pos]); continue;
2768 case kmax : pos--; tab[pos-1] = TMath::Max(tab[pos-1],tab[pos]); continue;
2769
2770 case klog : if (tab[pos-1] > 0) tab[pos-1] = TMath::Log(tab[pos-1]);
2771 else {tab[pos-1] = 0;} //{indetermination }
2772 continue;
2773 case kexp : { Double_t dexp = tab[pos-1];
2774 if (dexp < -700) {tab[pos-1] = 0; continue;}
2775 if (dexp > 700) {tab[pos-1] = TMath::Exp(700); continue;}
2776 tab[pos-1] = TMath::Exp(dexp); continue; }
2777 case klog10: if (tab[pos-1] > 0) tab[pos-1] = TMath::Log10(tab[pos-1]);
2778 else {tab[pos-1] = 0;} //{indetermination }
2779 continue;
2780
2781 case kpi : pos++; tab[pos-1] = TMath::ACos(-1); continue;
2782
2783 case kabs : tab[pos-1] = TMath::Abs(tab[pos-1]); continue;
2784 case ksign : if (tab[pos-1] < 0) tab[pos-1] = -1; else tab[pos-1] = 1; continue;
2785 case kint : tab[pos-1] = Double_t(Int_t(tab[pos-1])); continue;
2786
2787 case kSignInv: tab[pos-1] = -1 * tab[pos-1]; continue;
2788
2789 case krndm : pos++; tab[pos-1] = gRandom->Rndm(); continue;
2790
2791 case kAnd : pos--; if (tab[pos-1]!=0 && tab[pos]!=0) tab[pos-1]=1;
2792 else tab[pos-1]=0;
2793 continue;
2794 case kOr : pos--; if (tab[pos-1]!=0 || tab[pos]!=0) tab[pos-1]=1;
2795 else tab[pos-1]=0;
2796 continue;
2797 case kEqual: pos--; if (tab[pos-1] == tab[pos]) tab[pos-1]=1;
2798 else tab[pos-1]=0;
2799 continue;
2800 case kNotEqual : pos--; if (tab[pos-1] != tab[pos]) tab[pos-1]=1;
2801 else tab[pos-1]=0;
2802 continue;
2803 case kLess : pos--; if (tab[pos-1] < tab[pos]) tab[pos-1]=1;
2804 else tab[pos-1]=0;
2805 continue;
2806 case kGreater : pos--; if (tab[pos-1] > tab[pos]) tab[pos-1]=1;
2807 else tab[pos-1]=0;
2808 continue;
2809
2810 case kLessThan: pos--; if (tab[pos-1]<=tab[pos]) tab[pos-1]=1;
2811 else tab[pos-1]=0;
2812 continue;
2813 case kGreaterThan: pos--; if (tab[pos-1]>=tab[pos]) tab[pos-1]=1;
2814 else tab[pos-1]=0;
2815 continue;
2816 case kNot : if (tab[pos-1]!=0) tab[pos-1] = 0; else tab[pos-1] = 1;
2817 continue;
2818
2819 case kStringEqual : strpos -= 2; pos -=2 ; pos++;
2820 if (!strcmp(stringStack[strpos+1],stringStack[strpos])) tab[pos-1]=1;
2821 else tab[pos-1]=0;
2822 continue;
2823 case kStringNotEqual: strpos -= 2; pos -= 2; pos++;
2824 if (strcmp(stringStack[strpos+1],stringStack[strpos])) tab[pos-1]=1;
2825 else tab[pos-1]=0;
2826 continue;
2827
2828 case kBitAnd : pos--; tab[pos-1]= ((Int_t) tab[pos-1]) & ((Int_t) tab[pos]); continue;
2829 case kBitOr : pos--; tab[pos-1]= ((Int_t) tab[pos-1]) | ((Int_t) tab[pos]); continue;
2830 case kLeftShift : pos--; tab[pos-1]= ((Int_t) tab[pos-1]) <<((Int_t) tab[pos]); continue;
2831 case kRightShift: pos--; tab[pos-1]= ((Int_t) tab[pos-1]) >>((Int_t) tab[pos]); continue;
2832
2833 case kJump : i = (oper & kTFOperMask); continue;
2834 case kJumpIf : pos--; if (!tab[pos]) i = (oper & kTFOperMask); continue;
2835
2836 case kBoolOptimize: {
2837 // boolean operation optimizer
2838
2839 int param = (oper & kTFOperMask);
2840 Bool_t skip = kFALSE;
2841 int op = param % 10; // 1 is && , 2 is ||
2842
2843 if (op == 1 && (!tab[pos-1]) ) {
2844 // &&: skip the right part if the left part is already false
2845
2846 skip = kTRUE;
2847
2848 // Preserve the existing behavior (i.e. the result of a&&b is
2849 // either 0 or 1)
2850 tab[pos-1] = 0;
2851
2852 } else if (op == 2 && tab[pos-1] ) {
2853 // ||: skip the right part if the left part is already true
2854
2855 skip = kTRUE;
2856
2857 // Preserve the existing behavior (i.e. the result of a||b is
2858 // either 0 or 1)
2859 tab[pos-1] = 1;
2860
2861 }
2862
2863 if (skip) {
2864 int toskip = param / 10;
2865 i += toskip;
2866 }
2867 continue;
2868 }
2869
2870 }
2871
2872 switch(opcode) {
2873
2874 #define R__EXPO(var) \
2875 { \
2876 pos++; int param = (oper & kTFOperMask); \
2877 tab[pos-1] = TMath::Exp(params[param]+params[param+1]*x[var]); \
2878 continue; \
2879 }
2880 // case kexpo:
2881 case kxexpo: R__EXPO(0);
2882 case kyexpo: R__EXPO(1);
2883 case kzexpo: R__EXPO(2);
2884 case kxyexpo:{ pos++; int param = (oper & kTFOperMask);
2885 tab[pos-1] = TMath::Exp(params[param]+params[param+1]*x[0]+params[param+2]*x[1]);
2886 continue; }
2887
2888 #define R__GAUS(var) \
2889 { \
2890 pos++; int param = (oper & kTFOperMask); \
2891 tab[pos-1] = params[param]*TMath::Gaus(x[var],params[param+1],params[param+2],IsNormalized()); \
2892 continue; \
2893 }
2894
2895 // case kgaus:
2896 case kxgaus: R__GAUS(0);
2897 case kygaus: R__GAUS(1);
2898 case kzgaus: R__GAUS(2);
2899 case kxygaus: {pos++; int param = (oper & kTFOperMask);
2900 Double_t intermede1;
2901 if (params[param+2] == 0) {
2902 intermede1=1e10;
2903 } else {
2904 intermede1=Double_t((x[0]-params[param+1])/params[param+2]);
2905 }
2906 Double_t intermede2;
2907 if (params[param+4] == 0) {
2908 intermede2=1e10;
2909 } else {
2910 intermede2=Double_t((x[1]-params[param+3])/params[param+4]);
2911 }
2912 tab[pos-1] = params[param]*TMath::Exp(-0.5*(intermede1*intermede1+intermede2*intermede2));
2913 continue; }
2914
2915 #define R__LANDAU(var) \
2916 { \
2917 pos++; const int param = (oper & kTFOperMask); \
2918 tab[pos-1] = params[param]*TMath::Landau(x[var],params[param+1],params[param+2],IsNormalized()); \
2919 continue; \
2920 }
2921 // case klandau:
2922 case kxlandau: R__LANDAU(0);
2923 case kylandau: R__LANDAU(1);
2924 case kzlandau: R__LANDAU(2);
2925 case kxylandau: { pos++; int param = oper&0x7fffff /* ActionParams[i] */ ;
2926 Double_t intermede1=TMath::Landau(x[0], params[param+1], params[param+2],IsNormalized());
2927 Double_t intermede2=TMath::Landau(x[1], params[param+3], params[param+4],IsNormalized());
2928 tab[pos-1] = params[param]*intermede1*intermede2;
2929 continue;
2930 }
2931
2932 #define R__POLY(var) \
2933 { \
2934 pos++; int param = (oper & kTFOperMask); \
2935 tab[pos-1] = 0; Double_t intermede = 1; \
2936 Int_t inter = param/100; /* arrondit */ \
2937 Int_t int1= param-inter*100-1; /* aucune simplification ! (sic) */ \
2938 for (j=0 ;j<inter+1;j++) { \
2939 tab[pos-1] += intermede*params[j+int1]; \
2940 intermede *= x[var]; \
2941 } \
2942 continue; \
2943 }
2944 // case kpol:
2945 case kxpol: R__POLY(0);
2946 case kypol: R__POLY(1);
2947 case kzpol: R__POLY(2);
2948
2949 case kDefinedVariable : {
2950 if (!precalculated) {
2951 precalculated = 1;
2952 for(j=0;j<fNval;j++) param_calc[j]=DefinedValue(j);
2953 }
2954 pos++; tab[pos-1] = param_calc[(oper & kTFOperMask)];
2955 continue;
2956 }
2957
2958 case kDefinedString : {
2959 int param = (oper & kTFOperMask);
2960 if (!precalculated_str) {
2961 precalculated_str=1;
2962 for (j=0;j<fNstring;j++) string_calc[j]=DefinedString(j);
2963 }
2964 strpos++; stringStack[strpos-1] = string_calc[param];
2965 pos++; tab[pos-1] = 0;
2966 continue;
2967 }
2968
2969 case kFunctionCall: {
2970 // an external function call
2971
2972 int param = (oper & kTFOperMask);
2973 int fno = param / 1000;
2974 int nargs = param % 1000;
2975
2976 // Retrieve the function
2977 TMethodCall *method = (TMethodCall*)fFunctions.At(fno);
2978
2979 // Set the arguments
2980 method->ResetParam();
2981 if (nargs) {
2982 UInt_t argloc = pos-nargs;
2983 for(j=0;j<nargs;j++,argloc++,pos--) {
2984 method->SetParam(tab[argloc]);
2985 }
2986 }
2987 pos++;
2988 Double_t ret;
2989 method->Execute(ret);
2990 tab[pos-1] = ret; // check for the correct conversion!
2991
2992 continue;
2993 };
2994 }
2995 if (!TestBit(kOptimizationError)) {
2996 SetBit(kOptimizationError);
2997 Warning("EvalParOld","Found an unsupported opcode (%d)",oper >> kTFOperShift);
2998 }
2999 }
3000 Double_t result0 = tab[0];
3001 return result0;
3002
3003}
3004
3005////////////////////////////////////////////////////////////////////////////////
3006/// Reconstruct the formula expression from the internal TFormula member variables
3007///
3008/// This function uses the internal member variables of TFormula to
3009/// construct the mathematical expression associated with the TFormula
3010/// instance. This function can be used to get an expanded version of the
3011/// expression originally assigned to the TFormula instance, i.e. that
3012/// the string returned by GetExpFormula() doesn't depend on other
3013/// TFormula object names.
3014///
3015/// if option contains "p" the returned string will contain the formula
3016/// expression with symbolic parameters, eg [0] replaced by the actual value
3017/// of the parameter. Example:
3018/// if expression in formula is: "[0]*(x>-[1])+[2]*exp(-[3]*x)"
3019/// and parameters are 3.25,-4.01,4.44,-0.04, GetExpFormula("p") will return:
3020/// "(3.25*(x>+4.01))+(4.44*exp(+0.04*x))"
3021
3022TString TFormula::GetExpFormula(Option_t *option) const
3023{
3024 if (fNoper>0) {
3025 TString* tab=new TString[fNoper];
3026 Bool_t* ismulti=new Bool_t[fNoper];
3027 Int_t spos=0;
3028
3029 ismulti[0]=kFALSE;
3030 Int_t optype;
3031 Int_t j;
3032 Int_t ternaryend = -1;
3033 for(Int_t i=0;i<fNoper;i++){
3034 optype = GetAction(i);
3035
3036 if (ternaryend==i) {
3037 // The ? and : have been added to tab[spos-2]
3038 if(ismulti[spos-1]){
3039 tab[spos-2]=tab[spos-2]+"("+tab[spos-1]+")";
3040 } else {
3041 tab[spos-2]=tab[spos-2]+tab[spos-1];
3042 }
3043 spos--;
3044 // Do not call continue since we need to
3045 // do the rest of the loop.
3046 }
3047
3048 // Boolean optimization breakpoint
3049 if (optype==kBoolOptimize) { // -3) {
3050 continue;
3051 }
3052
3053 //Sign inversion
3054 if (optype==kSignInv) { // -1) {
3055 tab[spos-1]="-("+tab[spos-1]+")";
3056 // i++;
3057 continue;
3058 }
3059
3060 //Simple name (parameter,pol0,landau, etc)
3061 if (kexpo<=optype && optype<=kzpol) { // >=0) {
3062 tab[spos]=fExpr[i];
3063 ismulti[spos]=kFALSE;
3064 spos++;
3065 continue;
3066 }
3067
3068 //constants, variables x,y,z,t, pi
3069 if ((optype<=151 && optype>=140 && optype!=145) || (optype == 40)) {
3070 tab[spos]=fExpr[i];
3071 ismulti[spos]=kFALSE;
3072 spos++;
3073 continue;
3074 }
3075
3076 //Basic operators (+,-,*,/,==,^,etc)
3077 if(((optype>0 && optype<6) || optype==20 ||
3078 (((optype>59 && optype<69) || (optype >75 && optype<82)) && spos>=2))) {
3079 // if(optype==-20 && spos>=2){
3080 if(ismulti[spos-2]){
3081 tab[spos-2]="("+tab[spos-2]+")";
3082 }
3083 if(ismulti[spos-1]){
3084 tab[spos-2]+=fExpr[i]+("("+tab[spos-1]+")");
3085 }else{
3086 tab[spos-2]+=fExpr[i]+tab[spos-1];
3087 }
3088 ismulti[spos-2]=kTRUE;
3089 spos--;
3090 continue;
3091 }
3092 //Ternary condition
3093 if (optype==kJumpIf) {
3094 if(ismulti[spos-1]){
3095 tab[spos-1]="("+tab[spos-1]+")?";
3096 } else {
3097 tab[spos-1]=tab[spos-1]+"?";
3098 }
3099 continue;
3100 }
3101 if (optype==kJump) {
3102 if(ismulti[spos-1]){
3103 tab[spos-2]=tab[spos-2]+"("+tab[spos-1]+"):";
3104 } else {
3105 tab[spos-2]=tab[spos-2]+tab[spos-1]+":";
3106 }
3107 ternaryend = GetActionParam(i) + 1;
3108 spos--;
3109 continue;
3110 }
3111
3112 //Functions
3113 int offset = 0;
3114 TString funcname = fExpr[i];
3115 if((optype>9 && optype<16) ||
3116 (optype>20 && optype<23) ||
3117 (optype>29 && optype<34) ||
3118 (optype>40 && optype<44) ||
3119 (optype>69 && optype<76) ||
3120 (optype==145)) {
3121 //Functions with the format func(x)
3122 offset = -1;
3123 }
3124
3125 if((optype>15 && optype<20) ||
3126 (optype>22 && optype<26)) {
3127 //Functions with the format func(x,y)
3128 offset = -2;
3129 }
3130 if(optype==145) {
3131 int param = (fOper[i] & kTFOperMask);
3132 //int fno = param / 1000;
3133 int nargs = param % 1000;
3134 offset = -nargs;
3135 // The function name contains return type and parameters types we need
3136 // to trim them.
3137 int depth;
3138 for(j=0, depth=0;j<funcname.Length();++j) {
3139 switch (funcname[j]) {
3140 case '<':
3141 ++depth; break;
3142 case '>':
3143 --depth; break;
3144 case ' ':
3145 if (depth==0) {
3146 funcname.Remove(0,j+1);
3147 j = funcname.Length();
3148 break;
3149 }
3150 }
3151 }
3152 Ssiz_t ind = funcname.First('(');
3153 funcname.Remove(ind);
3154 }
3155 if (offset > 0) {
3156 Error("GetExpFormula","Internal error, number of argument found is %d",-offset);
3157 } else if (offset == 0) {
3158 tab[spos]=funcname+"()";
3159 ismulti[spos]=kFALSE;
3160 spos += 1;
3161 continue;
3162 } else if (offset<=0 && (spos+offset>=0)) {
3163 tab[spos+offset]=funcname+("("+tab[spos+offset]);
3164 for (j=offset+1; j<0; j++){
3165 tab[spos+offset]+=","+tab[spos+j];
3166 }
3167 tab[spos+offset]+=")";
3168 ismulti[spos+offset]=kFALSE;
3169 spos += offset+1;
3170 continue;
3171 }
3172 }
3173 if (ternaryend==fNoper) {
3174 // The ? and : have been added to tab[spos-2]
3175 if(ismulti[spos-1]){
3176 tab[spos-2]=tab[spos-2]+"("+tab[spos-1]+")";
3177 } else {
3178 tab[spos-2]=tab[spos-2]+tab[spos-1];
3179 }
3180 spos--;
3181 }
3182
3183 TString ret = "";
3184 if (spos > 0) ret = tab[spos-1];
3185 delete[] tab;
3186 delete[] ismulti;
3187
3188 //if option "p" is specified, return the real values of parameters instead of [0]
3189 TString opt = option;
3190 opt.ToLower();
3191 if (opt.Contains("p")) {
3192 char pb[13];
3193 char pbv[100];
3194 for (j=0;j<fNpar;j++) {
3195 snprintf(pb,sizeof(pb),"[%d]",j);
3196 snprintf(pbv,100,"%g",fParams[j]);
3197 ret.ReplaceAll(pb,pbv);
3198 }
3199 ret.ReplaceAll("--","+");
3200 ret.ReplaceAll("+-","-");
3201 }
3202 return ret;
3203 } else{
3204 TString ret="";
3205 return ret;
3206 }
3207}
3208
3209////////////////////////////////////////////////////////////////////////////////
3210/// Return linear part.
3211
3212const TObject* TFormula::GetLinearPart(Int_t i)
3213{
3214 if (!fLinearParts.IsEmpty())
3215 return fLinearParts.UncheckedAt(i);
3216 return 0;
3217}
3218
3219////////////////////////////////////////////////////////////////////////////////
3220/// Return value of parameter number ipar.
3221
3222Double_t TFormula::GetParameter(Int_t ipar) const
3223{
3224 if (ipar <0 || ipar >= fNpar) return 0;
3225 return fParams[ipar];
3226}
3227
3228////////////////////////////////////////////////////////////////////////////////
3229/// Return value of parameter named parName.
3230
3231Double_t TFormula::GetParameter(const char *parName) const
3232{
3233 const Double_t kNaN = 1e-300;
3234 Int_t index = GetParNumber(parName);
3235 if (index==-1) {
3236 Error("TFormula", "Parameter %s not found", parName);
3237 return kNaN;
3238 }
3239 return GetParameter(index);
3240}
3241
3242////////////////////////////////////////////////////////////////////////////////
3243/// Return name of one parameter.
3244
3245const char *TFormula::GetParName(Int_t ipar) const
3246{
3247 if (ipar <0 || ipar >= fNpar) return "";
3248 if (fNames[ipar].Length() > 0) return (const char*)fNames[ipar];
3249 return Form("p%d",ipar);
3250}
3251
3252////////////////////////////////////////////////////////////////////////////////
3253/// Return parameter number by name.
3254
3255Int_t TFormula::GetParNumber(const char *parName) const
3256{
3257 if (!parName)
3258 return -1;
3259
3260 for (Int_t i=0; i<fNpar; i++) {
3261 if (!strcmp(GetParName(i),parName)) return i;
3262 }
3263 return -1;
3264}
3265
3266////////////////////////////////////////////////////////////////////////////////
3267/// Return true if the expression at the index 'oper' has to be treated as a string
3268
3269Bool_t TFormula::IsString(Int_t oper) const
3270{
3271 return GetAction(oper) == kStringConst;
3272}
3273
3274////////////////////////////////////////////////////////////////////////////////
3275/// Dump this formula with its attributes.
3276
3277void TFormula::Print(Option_t *) const
3278{
3279 Int_t i;
3280 Printf(" %20s : %s Ndim= %d, Npar= %d, Noper= %d",GetName(),GetTitle(), fNdim,fNpar,fNoper);
3281 for (i=0;i<fNoper;i++) {
3282 Printf(" fExpr[%d] = %s action = %d action param = %d ",
3283 i,(const char*)fExpr[i],GetAction(i),GetActionParam(i));
3284 }
3285 //MI change
3286 //
3287 if (fNOperOptimized>0){
3288 Printf("Optimized expression");
3289 for (i=0;i<fNOperOptimized;i++) {
3290 Printf(" fExpr[%d] = %s\t\t action = %d action param = %d ",
3291 i,(const char*)fExprOptimized[i],GetActionOptimized(i),GetActionParamOptimized(i));
3292 }
3293 }
3294
3295 if (!fNames) return;
3296 if (!fParams) return;
3297 for (i=0;i<fNpar;i++) {
3298 Printf(" Par%3d %20s = %g",i,GetParName(i),fParams[i]);
3299 }
3300}
3301
3302////////////////////////////////////////////////////////////////////////////////
3303/// If the formula is for linear fitting, change the title to
3304/// normal and fill the LinearParts array
3305
3306void TFormula::ProcessLinear(TString &formula)
3307{
3308 TString formula2(formula);
3309 char repl[20];
3310 char *pch;
3311 Int_t nf, offset, replsize;
3312 //replace "++" with "+[i]*"
3313 pch= (char*)strstr(formula.Data(), "++");
3314 if (pch)
3315 formula.Insert(0, "[0]*(");
3316 pch= (char*)strstr(formula.Data(), "++");
3317 if (pch){
3318 //if there are "++", replaces them with +[i]*
3319 nf = 1;
3320 while (pch){
3321 snprintf(repl,20, ")+[%d]*(", nf);
3322 offset = pch-formula.Data();
3323 if (nf<10) replsize = 7;
3324 else if (nf<100) replsize = 8;
3325 else replsize = 9;
3326 formula.Replace(pch-formula.Data(), 2, repl, replsize);
3327 pch = (char*)strstr(formula.Data()+offset, "++");
3328 nf++;
3329 }
3330 formula.Append(')', 1);
3331 } else {
3332 //if there are no ++, create a new string with ++ instead of +[i]*
3333 formula2=formula2(4, formula2.Length()-4);
3334 pch= (char*)strchr(formula2.Data(), '[');
3335 snprintf(repl,20, "++");
3336 nf = 1;
3337 while (pch){
3338 offset = pch-formula2.Data()-1;
3339 if (nf<10) replsize = 5;
3340 else replsize = 6;
3341 formula2.Replace(pch-formula2.Data()-1, replsize, repl, 2);
3342 pch = (char*)strchr(formula2.Data()+offset, '[');
3343 nf++;
3344 }
3345 }
3346
3347 fLinearParts.Expand(nf);
3348 //break up the formula and fill the array of linear parts
3349 TString replaceformula;
3350 formula2 = formula2.ReplaceAll("++", 2, "|", 1);
3351 TObjArray *oa = formula2.Tokenize("|");
3352 TString replaceformula_name;
3353 for (Int_t i=0; i<nf; i++) {
3354 replaceformula = ((TObjString *)oa->UncheckedAt(i))->GetString();
3355 replaceformula_name = "f_linear_";
3356 replaceformula_name.Append(replaceformula);
3357 TFormula *f = new TFormula(replaceformula_name.Data(), replaceformula.Data());
3358 if (!f) {
3359 Error("TFormula", "f_linear not allocated");
3360 return;
3361 }
3362 {
3363 R__LOCKGUARD(gROOTMutex);
3364 gROOT->GetListOfFunctions()->Remove(f);
3365 }
3366 f->SetBit(kNotGlobal, 1);
3367 fLinearParts.Add(f);
3368 }
3369 oa->Delete();
3370}
3371
3372////////////////////////////////////////////////////////////////////////////////
3373/// Initialize parameter number ipar.
3374
3375void TFormula::SetParameter(const char *name, Double_t value)
3376{
3377 Int_t ipar = GetParNumber(name);
3378 if (ipar <0 || ipar >= fNpar) return;
3379 fParams[ipar] = value;
3380 Update();
3381}
3382
3383////////////////////////////////////////////////////////////////////////////////
3384/// Initialize parameter number ipar.
3385
3386void TFormula::SetParameter(Int_t ipar, Double_t value)
3387{
3388 if (ipar <0 || ipar >= fNpar) return;
3389 fParams[ipar] = value;
3390 Update();
3391}
3392
3393////////////////////////////////////////////////////////////////////////////////
3394/// Initialize array of all parameters.
3395/// See also the next function with the same name.
3396
3397void TFormula::SetParameters(const Double_t *params)
3398{
3399 for (Int_t i=0; i<fNpar;i++) {
3400 fParams[i] = params[i];
3401 }
3402 Update();
3403}
3404
3405////////////////////////////////////////////////////////////////////////////////
3406/// Initialize up to 11 parameters
3407/// All arguments except THE FIRST TWO are optional
3408/// In case of a function with only one parameter, call this function with p1=0.
3409/// Minimum two arguments are required to differentiate this function
3410/// from the SetParameters(cont Double_t *params)
3411
3412void TFormula::SetParameters(Double_t p0,Double_t p1,Double_t p2,Double_t p3,Double_t p4
3413 ,Double_t p5,Double_t p6,Double_t p7,Double_t p8,Double_t p9,Double_t p10)
3414{
3415 if (fNpar > 0) fParams[0] = p0;
3416 if (fNpar > 1) fParams[1] = p1;
3417 if (fNpar > 2) fParams[2] = p2;
3418 if (fNpar > 3) fParams[3] = p3;
3419 if (fNpar > 4) fParams[4] = p4;
3420 if (fNpar > 5) fParams[5] = p5;
3421 if (fNpar > 6) fParams[6] = p6;
3422 if (fNpar > 7) fParams[7] = p7;
3423 if (fNpar > 8) fParams[8] = p8;
3424 if (fNpar > 9) fParams[9] = p9;
3425 if (fNpar >10) fParams[10]= p10;
3426 Update();
3427}
3428
3429////////////////////////////////////////////////////////////////////////////////
3430/// Set name of parameter number ipar
3431
3432void TFormula::SetParName(Int_t ipar, const char *name)
3433{
3434 if (ipar <0 || ipar >= fNpar) return;
3435 fNames[ipar] = name;
3436}
3437
3438////////////////////////////////////////////////////////////////////////////////
3439/// Set up to 11 parameter names.
3440
3441void TFormula::SetParNames(const char*name0,const char*name1,const char*name2,const char*name3,const char*name4,
3442 const char*name5,const char*name6,const char*name7,const char*name8,const char*name9,const char*name10)
3443{
3444 if (fNpar > 0) fNames[0] = name0;
3445 if (fNpar > 1) fNames[1] = name1;
3446 if (fNpar > 2) fNames[2] = name2;
3447 if (fNpar > 3) fNames[3] = name3;
3448 if (fNpar > 4) fNames[4] = name4;
3449 if (fNpar > 5) fNames[5] = name5;
3450 if (fNpar > 6) fNames[6] = name6;
3451 if (fNpar > 7) fNames[7] = name7;
3452 if (fNpar > 8) fNames[8] = name8;
3453 if (fNpar > 9) fNames[9] = name9;
3454 if (fNpar >10) fNames[10]= name10;
3455}
3456
3457////////////////////////////////////////////////////////////////////////////////
3458/// Stream a class object.
3459
3460void TFormula::Streamer(TBuffer &b, const TClass *onfile_class)
3461{
3462 if (b.IsReading()) {
3463 UInt_t R__s, R__c;
3464 Version_t v = b.ReadVersion(&R__s, &R__c);
3465 if (v==6) {
3466 Error("Streamer","version 6 is not supported");
3467 return;
3468 }
3469 Streamer(b, v, R__s, R__c, onfile_class);
3470
3471 } else {
3472 b.WriteClassBuffer(TFormula::Class(),this);
3473 }
3474}
3475
3476////////////////////////////////////////////////////////////////////////////////
3477/// Stream a class object.
3478
3479void TFormula::Streamer(TBuffer &b)
3480{
3481 if (b.IsReading()) {
3482 UInt_t R__s, R__c;
3483 Version_t v = b.ReadVersion(&R__s, &R__c);
3484 if (v==6) {
3485 Error("Streamer","version 6 is not supported");
3486 return;
3487 }
3488 Streamer(b, v, R__s, R__c, nullptr);
3489
3490 } else {
3491 b.WriteClassBuffer(TFormula::Class(),this);
3492 }
3493}
3494
3495////////////////////////////////////////////////////////////////////////////////
3496/// specialized streamer function being able to read old TF1 versions as TF1Old in memory
3497
3498void TFormula::Streamer(TBuffer &b, Int_t v, UInt_t R__s, UInt_t R__c, const TClass *onfile_class)
3499{
3500
3501 //printf("Reading TFormula - version %d \n",v);
3502 if (v > 3 ) {
3503 b.ReadClassBuffer(TFormula::Class(), this, v, R__s, R__c, onfile_class);
3504 if (!TestBit(kNotGlobal)) {
3505 R__LOCKGUARD(gROOTMutex);
3506 gROOT->GetListOfFunctions()->Add(this);
3507 }
3508
3509 // We need to reinstate (if possible) the TMethodCall.
3510 if (fFunctions.GetLast()>=0) {
3511 // Compiles will reset the parameter values so we need
3512 // to temporarily keep them
3513 Double_t *param = fParams;
3514 TString *names = fNames;
3515 Int_t npar = fNpar;
3516 fParams = 0;
3517 fNames = 0;
3518 if (Compile()) {
3519 Error("Streamer","error compiling formula");
3520 return;
3521 }
3522 for (Int_t i = 0; i<npar && i<fNpar; ++i) fParams[i] = param[i];
3523 delete [] param;
3524 delete [] fNames;
3525 fNames = names;
3526 } else if (v<6) {
3527 Convert(v);
3528 }
3529 Optimize();
3530 return;
3531 }
3532 // version smaller or equal to 3
3533 // process old versions before automatic schema evolution
3534 TNamed::Streamer(b);
3535 b >> fNdim;
3536 b >> fNumber;
3537 if (v > 1) b >> fNval;
3538 if (v > 2) b >> fNstring;
3539 fNpar = b.ReadArray(fParams);
3540 fOper = new Int_t[gMAXOP];
3541 fNoper = b.ReadArray(fOper);
3542 fNconst = b.ReadArray(fConst);
3543 if (fNoper) {
3544 fExpr = new TString[fNoper];
3545 }
3546 if (fNpar) {
3547 fNames = new TString[fNpar];
3548 }
3549 Int_t i;
3550 for (i=0;i<fNoper;i++) fExpr[i].Streamer(b);
3551 for (i=0;i<fNpar;i++) fNames[i].Streamer(b);
3552 {
3553 R__LOCKGUARD(gROOTMutex);
3554 if (gROOT->GetListOfFunctions()->FindObject(GetName())) return;
3555 gROOT->GetListOfFunctions()->Add(this);
3556 }
3557 b.CheckByteCount(R__s, R__c, TFormula::IsA());
3558
3559 Convert(v);
3560 // end of old versions
3561
3562}
3563
3564void TFormula::Convert(UInt_t /* fromVersion */)
3565{
3566 // Convert the fOper of a TFormula version fromVersion to the current in memory version
3567
3568 enum {
3569 kOldexpo = 1000,
3570 kOldgaus = 2000,
3571 kOldlandau = 4000,
3572 kOldxylandau = 4500,
3573 kOldConstants = 50000,
3574 kOldStrings = 80000,
3575 kOldVariable = 100000,
3576 kOldTreeString = 105000,
3577 kOldFormulaVar = 110000,
3578 kOldBoolOptimize = 120000,
3579 kOldFunctionCall = 200000
3580 };
3581 int i,j;
3582
3583 for (i=0,j=0; i<fNoper; ++i,++j) {
3584 Int_t action = fOper[i];
3585 Int_t newActionCode = 0;
3586 Int_t newActionParam = 0;
3587
3588 if ( action == 0) {
3589 // Sign Inversion
3590
3591 newActionCode = kSignInv;
3592
3593 Float_t aresult = 99.99;
3594 sscanf((const char*)fExpr[i],"%g",&aresult);
3595 R__ASSERT((aresult+1)<0.001);
3596
3597 ++i; // skip the implied multiplication.
3598
3599 // For consistency and for Optimize to work correctly
3600 // we need to remove the "-1" string in fExpr
3601 for (int z=i; z<fNoper; ++z) {
3602 fExpr[z-1] = fExpr[z];
3603 }
3604
3605 } else if ( action < 100 ) {
3606 // basic operators and mathematical library
3607
3608 newActionCode = action;
3609
3610 } else if (action >= kOldFunctionCall) {
3611 // Function call
3612
3613 newActionCode = kFunctionCall;
3614 newActionParam = action-kOldFunctionCall;
3615
3616 } else if (action >= kOldBoolOptimize) {
3617 // boolean operation optimizer
3618
3619 newActionCode = kBoolOptimize;
3620 newActionParam = action-kOldBoolOptimize;
3621
3622 } else if (action >= kOldFormulaVar) {
3623 // a variable
3624
3625 newActionCode = kVariable;
3626 newActionParam = action-kOldFormulaVar;
3627
3628 } else if (action >= kOldTreeString) {
3629 // a tree string
3630
3631 newActionCode = kDefinedString;
3632 newActionParam = action-kOldTreeString;
3633
3634 } else if (action >= kOldVariable) {
3635 // a tree variable
3636
3637 newActionCode = kDefinedVariable;
3638 newActionParam = action-kOldVariable;
3639
3640 } else if (action == kOldStrings) {
3641 // String
3642
3643 newActionCode = kStringConst;
3644
3645 } else if (action >= kOldConstants) {
3646 // numerical value
3647
3648 newActionCode = kConstant;
3649 newActionParam = action-kOldConstants;
3650
3651 } else if (action > 10000 && action < kOldConstants) {
3652 // Polynomial
3653
3654 int var = action/10000; //arrondit
3655 newActionCode = kpol + (var-1);
3656 newActionParam = action - var*10000;
3657
3658 } else if (action >= 4600) {
3659
3660 Error("Convert","Unsupported value %d",action);
3661
3662 } else if (action > kOldxylandau) {
3663 // xylandau
3664
3665 newActionCode = kxylandau;
3666 newActionParam = action - (kOldxylandau+1);
3667
3668 } else if (action > kOldlandau) {
3669 // landau, xlandau, ylandau or zlandau
3670
3671 newActionCode = klandau;
3672 int var = action/100-40;
3673 if (var) newActionCode += var;
3674 newActionParam = action - var*100 - (kOldlandau+1);
3675
3676 } else if (action > 2500 && action < 2600) {
3677 // xygaus
3678
3679 newActionCode = kxygaus;
3680 newActionParam = action-2501;
3681
3682 } else if (action > 2000 && action < 2500) {
3683 // gaus, xgaus, ygaus or zgaus
3684
3685 newActionCode = kgaus;
3686 int var = action/100-20;
3687 if (var) newActionCode += var;
3688 newActionParam = action - var*100 - (kOldgaus+1);
3689
3690 } else if (action > 1500 && action < 1600) {
3691 // xyexpo
3692
3693 newActionCode = kxyexpo;
3694 newActionParam = action-1501;
3695
3696 } else if (action > 1000 && action < 1500) {
3697 // expo or xexpo or yexpo or zexpo
3698
3699 newActionCode = kexpo;
3700 int var = action/100-10;
3701 if (var) newActionCode += var;
3702 newActionParam = action - var*100 - (kOldexpo+1);
3703
3704 } else if (action > 100 && action < 200) {
3705 // Parameter substitution
3706
3707 newActionCode = kParameter;
3708 newActionParam = action - 101;
3709 }
3710
3711 SetAction( j, newActionCode, newActionParam );
3712
3713 }
3714 if (i!=j) {
3715 fNoper -= (i-j);
3716 }
3717
3718}
3719
3720////////////////////////////////////////////////////////////////////////////////
3721/// TOper offset - helper class for TFormula*
3722/// specify type of operand
3723/// fTypeX = kVariable
3724/// = kParameter
3725/// = kConstant
3726/// fOffestX = offset in corresponding array
3727
3728TOperOffset::TOperOffset()
3729{
3730 fType0=0;
3731 fType1=0;
3732 fType2=0;
3733 fType3=0;
3734 fOffset0=0;
3735 fOffset1=0;
3736 fOffset2=0;
3737 fOffset3=0;
3738 fOldAction=0;
3739 fToJump=0;
3740}
3741
3742////////////////////////////////////////////////////////////////////////////////
3743/// MakePrimitive
3744/// find TFormulaPrimitive replacement for some operands
3745
3746void TFormula::MakePrimitive(const char *expr, Int_t pos)
3747{
3748 TString cbase(expr);
3749 cbase.ReplaceAll("Double_t ","");
3750 int paran = cbase.First("(");
3751 // int nargs = 0;
3752 if (paran>0) {
3753 //nargs = 1;
3754 cbase[paran]=0;
3755 }
3756
3757 if (cbase=="<") cbase="XlY";
3758 if (cbase=="<=") cbase="XleY";
3759 if (cbase==">") cbase="XgY";
3760 if (cbase==">=") cbase="XgeY";
3761 if (cbase=="==" && GetActionOptimized(pos)!=kStringEqual) cbase="XeY";
3762 if (cbase=="!=" && GetActionOptimized(pos)!=kStringNotEqual) cbase="XneY";
3763
3764 ROOT::v5::TFormulaPrimitive *prim = ROOT::v5::TFormulaPrimitive::FindFormula(cbase ,paran>0 ? cbase.Data() + paran + 1 : (const char*)0);
3765 if (prim) {
3766 fPredefined[pos] = prim;
3767 if (prim->fType==10) {
3768 SetActionOptimized(pos, kFD1);
3769 }
3770 if (prim->fType==110) {
3771 SetActionOptimized(pos, kFD2);
3772 }
3773 if (prim->fType==1110) {
3774 SetActionOptimized(pos, kFD3);
3775 }
3776 if (prim->fType==-1) {
3777 SetActionOptimized(pos, kFDM);
3778 }
3779 if (prim->fType==0){
3780 SetActionOptimized(pos,kConstant,fNconst);
3781 fConst[fNconst] = prim->Eval(0);
3782 fNconst++;
3783 }
3784 return;
3785 }
3786}
3787
3788////////////////////////////////////////////////////////////////////////////////
3789/// MI include
3790///
3791/// Optimize formula
3792/// - Minimize the number of operands
3793/// 1. several operands are glued together
3794/// 2. some primitive functions glued together - exemp. (x+y) => PlusXY(x,y)
3795/// 3. maximize number of standard calls minimizing number of jumps in Eval cases
3796/// 4. variables, parameters and constants are mapped - using fOperOfssets0
3797/// Eval procedure use direct acces to data (only one corresponding case statement in eval procedure)
3798/// ~~~ {.cpp}
3799/// pdata[operand={Var,Par,Const}][offset]
3800/// pdata[fOperOffsets0[i]][fOperOffset1[i+1]]
3801/// ~~~
3802/// - The fastest evaluation function is chosen at the end
3803/// 1. fOptimal := pointer to the fastest function for given evaluation string
3804/// ~~~ {.cpp}
3805/// switch(GetActionOptimized(0)){
3806/// case kData : {fOptimal= (TFormulaPrimitive::TFuncG)&TFormula::EvalPrimitive0; break;}
3807/// case kUnary : {fOptimal= (TFormulaPrimitive::TFuncG)&TFormula::EvalPrimitive1; break;}
3808/// case kBinary : {fOptimal= (TFormulaPrimitive::TFuncG)&TFormula::EvalPrimitive2; break;}
3809/// case kThree : {fOptimal= (TFormulaPrimitive::TFuncG)&TFormula::EvalPrimitive3; break;}
3810/// case kFDM : {fOptimal= (TFormulaPrimitive::TFuncG)&TFormula::EvalPrimitive4; break;}
3811/// }
3812/// ~~~
3813/// 2. ex.
3814/// - fOptimal = ::EvalPrimitive0 - if it return only variable, constant or parameter
3815/// - = ::EvalParameter1 - if only one unary operation
3816/// - = ::EvalPrimitive2 - if only one binary operation
3817
3818void TFormula::Optimize()
3819{
3820 //
3821 // Initialize data members
3822 //
3823
3824 Int_t i;
3825
3826 if (fPredefined) { delete [] fPredefined; fPredefined = 0;}
3827 if (fOperOffset) { delete [] fOperOffset; fOperOffset = 0;}
3828 if (fExprOptimized) { delete [] fExprOptimized; fExprOptimized = 0;}
3829 if (fOperOptimized) { delete [] fOperOptimized; fOperOptimized = 0;}
3830
3831 fExprOptimized = new TString[fNoper];
3832 fOperOptimized = new Int_t[fNoper];
3833 fPredefined = new ROOT::v5::TFormulaPrimitive*[fNoper];
3834 fOperOffset = new TOperOffset[fNoper];
3835 for (i=0; i<fNoper; i++) {
3836 fExprOptimized[i] = fExpr[i] ;
3837 fOperOptimized[i] = fOper[i];
3838 fPredefined[i]= 0;
3839 }
3840
3841 //
3842 //Make primitives
3843 //
3844 for (i=0;i<fNoper;i++){
3845 if (fExprOptimized[i].Data()) {
3846 MakePrimitive(fExprOptimized[i].Data(), i);
3847 }
3848 }
3849 //
3850 Int_t maxfound = fNoper+1;
3851 Int_t *offset = new Int_t[maxfound*16];
3852 Int_t *optimized = new Int_t[maxfound];
3853 //
3854 //
3855 ROOT::v5::TFormulaPrimitive* primitive[10];
3856 primitive[0] = ROOT::v5::TFormulaPrimitive::FindFormula("PlusXY");
3857 primitive[1] = ROOT::v5::TFormulaPrimitive::FindFormula("MinusXY");
3858 primitive[2] = ROOT::v5::TFormulaPrimitive::FindFormula("MultXY");
3859 primitive[3] = ROOT::v5::TFormulaPrimitive::FindFormula("DivXY");
3860 primitive[4] = ROOT::v5::TFormulaPrimitive::FindFormula("XpYpZ");
3861 primitive[5] = ROOT::v5::TFormulaPrimitive::FindFormula("XxYxZ");
3862 primitive[6] = ROOT::v5::TFormulaPrimitive::FindFormula("XxYpZ");
3863 primitive[7] = ROOT::v5::TFormulaPrimitive::FindFormula("XpYxZ");
3864 primitive[8] = ROOT::v5::TFormulaPrimitive::FindFormula("Pow2");
3865 primitive[9] = ROOT::v5::TFormulaPrimitive::FindFormula("Pow3");
3866 //
3867 // set data pointers
3868 //
3869 for (i=0;i<fNoper;i++) optimized[i]=0;
3870 //
3871 for (i=0;i<fNoper;i++){
3872 Int_t actionparam = GetActionParamOptimized(i);
3873 Int_t action = GetActionOptimized(i);
3874
3875 if (action==kBoolOptimize){
3876 //
3877 // optimize booleans
3878 //
3879 fOperOffset[i].fType1 = actionparam/10; // operands to skip
3880 fOperOffset[i].fOffset0 = actionparam%10; // 1 is && , 2 is || - operand
3881 fOperOffset[i].fToJump = i+fOperOffset[i].fType1; // where we should jump
3882 continue;
3883 }
3884 if (action==kJump || action==kJumpIf) {
3885 // Ternary conditional operator
3886 fOperOffset[i].fType1 = action;
3887 fOperOffset[i].fToJump = actionparam;
3888 }
3889 //
3890 if (action==kConstant&&i<fNoper-2){
3891 //
3892 // get offsets for kFDM operands
3893 //
3894 if (GetActionOptimized(i+1)==kConstant && GetActionOptimized(i+2)==kFDM){
3895 optimized[i]=1;
3896 optimized[i+1]=1;
3897 i+=2;
3898 fOperOffset[i].fType0=actionparam;
3899 fOperOffset[i].fOffset0=GetActionParamOptimized(i-1);
3900 Int_t offset2 = int(fConst[fOperOffset[i].fOffset0]+0.4);
3901 fOperOffset[i].fOffset0=offset2;
3902 Int_t nparmax = offset2+fPredefined[i]->fNParameters;
3903 if (nparmax>fNpar){ // increase expected number of parameters
3904 fNpar=nparmax;
3905 }
3906 continue;
3907 }
3908 }
3909 switch(action){
3910 case kVariable : {action=kData; fOperOffset[i].fType0=0; break;}
3911 case kParameter: {action=kData; fOperOffset[i].fType0=1; break;}
3912 case kConstant : {action=kData; fOperOffset[i].fType0=2; break;}
3913 }
3914 //
3915 fOperOffset[i].fOffset0 = GetActionParamOptimized(i);
3916 SetActionOptimized(i,action, actionparam); //set common data option
3917 }
3918 //
3919 //
3920 fNOperOptimized = fNoper;
3921 //
3922 for (i=0; i<fNoper; ++i)
3923 {
3924 //
3925 if (!(GetActionOptimized(i)== kData)) continue;
3926 offset[0] = fOperOffset[i].fType0; //
3927 offset[1] = fOperOptimized[i] & kTFOperMask; // offset
3928
3929 if ((i+1) >= fNoper) continue;
3930
3931 if (GetActionOptimized(i+1)==kFD1){
3932 optimized[i] = 1; // to be optimized
3933 i++;
3934 fOperOffset[i].fType0 = offset[0];
3935 fOperOffset[i].fOffset0 = offset[1];
3936 SetActionOptimized(i ,kUnary);
3937 continue;
3938 }
3939 if (GetActionOptimized(i+1)==kAdd){
3940 optimized[i] = 1; // to be optimized
3941 i++;
3942 fOperOffset[i].fType0 = offset[0];
3943 fOperOffset[i].fOffset0 = offset[1];
3944 SetActionOptimized(i ,kPlusD);
3945 continue;
3946 }
3947 if (GetActionOptimized(i+1)==kMultiply){
3948 optimized[i] = 1; // to be optimized
3949 i++;
3950 fOperOffset[i].fType0 = offset[0];
3951 fOperOffset[i].fOffset0 = offset[1];
3952 SetActionOptimized(i,kMultD);
3953 continue;
3954 }
3955
3956 if ((i+2) >= fNoper) continue;
3957
3958 //
3959 //Binary operators
3960 if (!(GetActionOptimized(i+1)== kData)) continue;
3961 offset[2] = fOperOffset[i+1].fType0;
3962 offset[3] = fOperOptimized[i+1] & kTFOperMask; // offset
3963 //
3964 if (GetActionOptimized(i+2)==kFD2 || GetActionOptimized(i+2)==kAdd ||GetActionOptimized(i+2)==kSubstract||
3965 GetActionOptimized(i+2)==kMultiply || GetActionOptimized(i+2)==kDivide){
3966
3967 optimized[i] = 1; // to be optimized
3968 optimized[i+1] = 1; // to be optimized
3969 i+=2;
3970 //
3971 fOperOffset[i].fType0 = offset[0];
3972 fOperOffset[i].fOffset0 = offset[1];
3973 fOperOffset[i].fType1 = offset[2];
3974 fOperOffset[i].fOffset1 = offset[3];
3975 fOperOffset[i].fType2 = GetActionOptimized(i); //remember old action
3976 if (GetActionOptimized(i)==kAdd) {fPredefined[i] = primitive[0];}
3977 if (GetActionOptimized(i)==kSubstract) {fPredefined[i] = primitive[1];}
3978 if (GetActionOptimized(i)==kMultiply) {
3979 fPredefined[i]=primitive[2];
3980 if (offset[0]==offset[2]&&offset[1]==offset[3]) {
3981 fPredefined[i] = primitive[8];
3982 SetActionOptimized(i,kUnary);
3983 continue;
3984 }
3985 }
3986 if (GetActionOptimized(i)==kDivide) {
3987 fPredefined[i] = primitive[3];
3988 }
3989 SetActionOptimized(i,kBinary);
3990 continue;
3991 }
3992
3993 if ((i+3) >= fNoper) continue;
3994
3995 //
3996 //operator 3
3997 //
3998 if (!(GetActionOptimized(i+2)== kData)) continue;
3999 offset[4] = fOperOffset[i+2].fType0;
4000 offset[5] = fOperOptimized[i+2] & kTFOperMask; // offset
4001 //
4002 if (GetActionOptimized(i+3)==kFD3|| ( (GetActionOptimized(i+3)==kAdd||GetActionOptimized(i+3)==kMultiply) &&
4003 (GetActionOptimized(i+4)==kAdd||GetActionOptimized(i+4)==kMultiply) ) ){
4004 optimized[i+0] = 1; // to be optimized
4005 optimized[i+1] = 1; // to be optimized
4006 optimized[i+2] = 1; // to be optimized
4007 i+=3;
4008 //
4009 fOperOffset[i].fType0 = offset[0];
4010 fOperOffset[i].fOffset0 = offset[1];
4011 fOperOffset[i].fType1 = offset[2];
4012 fOperOffset[i].fOffset1 = offset[3];
4013 fOperOffset[i].fType2 = offset[4];
4014 fOperOffset[i].fOffset2 = offset[5];
4015 //
4016 fOperOffset[i].fOldAction = GetActionOptimized(i); //remember old action
4017 if (GetActionOptimized(i)==kFD3) {
4018 SetActionOptimized(i,kThree);
4019 continue;
4020 }
4021 Int_t action=0;
4022 Int_t action2=kThree;
4023 if (GetActionOptimized(i)==kAdd&&GetActionOptimized(i+1)==kAdd) action=4;
4024 if (GetActionOptimized(i)==kMultiply&&GetActionOptimized(i+1)==kMultiply) {
4025 action=5;
4026 if (offset[0]==offset[2]&&offset[1]==offset[3]&&offset[0]==offset[4]&&offset[1]==offset[5]){
4027 fPredefined[i]=primitive[9];
4028 action2=kUnary;
4029 action =9;
4030 }
4031 }
4032 if (GetActionOptimized(i)==kAdd&&GetActionOptimized(i+1)==kMultiply) action=6;
4033 if (GetActionOptimized(i)==kMultiply&&GetActionOptimized(i+1)==kAdd) action=7;
4034 //
4035 optimized[i]=1;
4036 i++;
4037 fOperOffset[i].fType0 = offset[0];
4038 fOperOffset[i].fOffset0 = offset[1];
4039 fOperOffset[i].fType1 = offset[2];
4040 fOperOffset[i].fOffset1 = offset[3];
4041 fOperOffset[i].fType2 = offset[4];
4042 fOperOffset[i].fOffset2 = offset[5];
4043 fPredefined[i]=primitive[action];
4044 SetActionOptimized(i,action2);
4045 continue;
4046 }
4047 }
4048 //
4049 //
4050 Int_t operO=0;
4051 TString expr="";
4052 Int_t *map0 = new Int_t[maxfound]; //remapping of the operands
4053 Int_t *map1 = new Int_t[maxfound]; //remapping of the operands
4054 for (i=0;i<fNoper;i++){
4055 map0[i] = operO;
4056 map1[operO] = i;
4057 fOperOptimized[operO] = fOperOptimized[i];
4058 fPredefined[operO] = fPredefined[i];
4059 fOperOffset[operO] = fOperOffset[i];
4060 expr += fExprOptimized[i];
4061 if (optimized[i]==0){
4062 fExprOptimized[operO] = expr;
4063 expr = "";
4064 operO++;
4065 }else{
4066 expr += ",";
4067 }
4068 }
4069 //
4070 // Recalculate long jump for Boolean optimize
4071 //
4072 for (i=0; i<fNOperOptimized; i++){
4073 Int_t optaction = GetActionOptimized(i);
4074 if (optaction==kBoolOptimize){
4075 Int_t oldpos = fOperOffset[i].fToJump;
4076 Int_t newpos = oldpos==fNoper ? fNOperOptimized : map0[oldpos];
4077 fOperOffset[i].fToJump = newpos; // new position to jump
4078 Int_t actionop = GetActionParamOptimized(i) % 10;
4079 switch (actionop) {
4080 case 1: SetActionOptimized(i,kBoolOptimizeAnd,newpos); break;
4081 case 2: SetActionOptimized(i,kBoolOptimizeOr,newpos); break;
4082 }
4083 } else if (optaction==kJump || optaction==kJumpIf) {
4084 Int_t oldpos = fOperOffset[i].fToJump;
4085 Int_t newpos = oldpos==fNoper ? fNOperOptimized : map0[oldpos];
4086 fOperOffset[i].fToJump = newpos; // new position to jump
4087 SetActionOptimized(i,optaction,newpos);
4088 }
4089 }
4090
4091
4092 fNOperOptimized = operO;
4093 //
4094 fOptimal= (ROOT::v5::TFormulaPrimitive::TFuncG)&TFormula::EvalParFast;
4095 if (fNOperOptimized==1) {
4096 switch(GetActionOptimized(0)){
4097 case kData : {fOptimal= (ROOT::v5::TFormulaPrimitive::TFuncG)&TFormula::EvalPrimitive0; break;}
4098 case kUnary : {fOptimal= (ROOT::v5::TFormulaPrimitive::TFuncG)&TFormula::EvalPrimitive1; break;}
4099 case kBinary : {fOptimal= (ROOT::v5::TFormulaPrimitive::TFuncG)&TFormula::EvalPrimitive2; break;}
4100 case kThree : {fOptimal= (ROOT::v5::TFormulaPrimitive::TFuncG)&TFormula::EvalPrimitive3; break;}
4101 case kFDM : {fOptimal= (ROOT::v5::TFormulaPrimitive::TFuncG)&TFormula::EvalPrimitive4; break;}
4102 }
4103 }
4104
4105 delete [] map1;
4106 delete [] map0;
4107 delete [] offset;
4108 delete [] optimized;
4109}
4110
4111////////////////////////////////////////////////////////////////////////////////
4112/// Evaluate primitive formula
4113
4114Double_t TFormula::EvalPrimitive(const Double_t *x, const Double_t *params)
4115{
4116 const Double_t *pdata[3] = {x,(params!=0)?params:fParams, fConst};
4117 Double_t result = pdata[fOperOffset->fType0][fOperOffset->fOffset0];
4118 switch((fOperOptimized[0] >> kTFOperShift)) {
4119 case kData : return result;
4120 case kUnary : return (fPredefined[0]->fFunc10)(pdata[fOperOffset->fType0][fOperOffset->fOffset0]);
4121 case kBinary :return (fPredefined[0]->fFunc110)(result,
4122 pdata[fOperOffset->fType1][fOperOffset->fOffset1]);
4123
4124 case kThree :return (fPredefined[0]->fFunc1110)(result, pdata[fOperOffset->fType1][fOperOffset->fOffset1],
4125 pdata[fOperOffset->fType2][fOperOffset->fOffset2]);
4126 case kFDM : return (fPredefined[0]->fFuncG)((Double_t*)&x[fOperOffset->fType0],
4127 (Double_t*)&params[fOperOffset->fOffset0]);
4128 }
4129 return 0;
4130}
4131
4132////////////////////////////////////////////////////////////////////////////////
4133/// Evaluate primitive formula
4134
4135Double_t TFormula::EvalPrimitive0(const Double_t *x, const Double_t *params)
4136{
4137 const Double_t *pdata[3] = {x,(params!=0)?params:fParams, fConst};
4138 return pdata[fOperOffset->fType0][fOperOffset->fOffset0];
4139}
4140
4141////////////////////////////////////////////////////////////////////////////////
4142/// Evaluate primitive formula
4143
4144Double_t TFormula::EvalPrimitive1(const Double_t *x, const Double_t *params)
4145{
4146 const Double_t *pdata[3] = {x,(params!=0)?params:fParams, fConst};
4147 return (fPredefined[0]->fFunc10)(pdata[fOperOffset->fType0][fOperOffset->fOffset0]);
4148}
4149
4150////////////////////////////////////////////////////////////////////////////////
4151/// Evaluate primitive formula
4152
4153Double_t TFormula::EvalPrimitive2(const Double_t *x, const Double_t *params)
4154{
4155 const Double_t *pdata[3] = {x,(params!=0)?params:fParams, fConst};
4156 return (fPredefined[0]->fFunc110)(pdata[fOperOffset->fType0][fOperOffset->fOffset0],
4157 pdata[fOperOffset->fType1][fOperOffset->fOffset1]);
4158}
4159
4160////////////////////////////////////////////////////////////////////////////////
4161/// Evaluate primitive formula
4162
4163Double_t TFormula::EvalPrimitive3(const Double_t *x, const Double_t *params)
4164{
4165 const Double_t *pdata[3] = {x,(params!=0)?params:fParams, fConst};
4166 return (fPredefined[0]->fFunc1110)(pdata[fOperOffset->fType0][fOperOffset->fOffset0], pdata[fOperOffset->fType1][fOperOffset->fOffset1],
4167 pdata[fOperOffset->fType2][fOperOffset->fOffset2]);
4168}
4169
4170////////////////////////////////////////////////////////////////////////////////
4171/// Evaluate primitive formula
4172
4173Double_t TFormula::EvalPrimitive4(const Double_t *x, const Double_t *params)
4174{
4175 const Double_t *par = (params!=0)?params:fParams;
4176 return (fPredefined[0]->fFuncG)((Double_t*)&x[fOperOffset->fType0],
4177 (Double_t*)&par[fOperOffset->fOffset0]);
4178}
4179
4180////////////////////////////////////////////////////////////////////////////////
4181/// Evaluate this formula.
4182///
4183/// The current value of variables x,y,z,t is passed through the pointer x.
4184/// The parameters used will be the ones in the array params if params is given
4185/// otherwise parameters will be taken from the stored data members fParams
4186///
4187/// \image html TFormula_eval.png
4188
4189Double_t TFormula::EvalParFast(const Double_t *x, const Double_t *uparams)
4190{
4191 const Double_t *pdata[3] = {x,(uparams!=0)?uparams:fParams, fConst};
4192 //
4193 Int_t i,j;
4194 Double_t tab[kMAXFOUND] = {0};
4195 const char *stringStack[gMAXSTRINGFOUND] = {0};
4196 Double_t param_calc[kMAXFOUND];
4197 char *string_calc[gMAXSTRINGFOUND] = {0};
4198 Int_t precalculated = 0;
4199 Int_t precalculated_str = 0;
4200
4201 Double_t *params;
4202
4203 if (uparams) {
4204 //for (j=0;j<fNpar;j++) fParams[j] = params[j];
4205 params = const_cast<Double_t*>(uparams);
4206 } else {
4207 params = fParams;
4208 }
4209
4210 //if (params) {
4211 // for (j=0;j<fNpar;j++) fParams[j] = params[j];
4212 //}
4213 UInt_t pos = 0;
4214 UInt_t strpos = 0;
4215 // for (i=0; i<fNoper; ++i) {
4216 for (i=0; i<fNOperOptimized; ++i) {
4217 //
4218 const int oper = fOperOptimized[i];
4219 const int opcode = oper >> kTFOperShift;
4220
4221 switch(opcode) { // FREQUENTLY USED OPERATION
4222 case kData : tab[pos] = pdata[fOperOffset[i].fType0][fOperOffset[i].fOffset0]; pos++;continue;
4223 case kPlusD : tab[pos-1]+= pdata[fOperOffset[i].fType0][fOperOffset[i].fOffset0]; continue;
4224 case kMultD : tab[pos-1]*= pdata[fOperOffset[i].fType0][fOperOffset[i].fOffset0]; continue;
4225 case kAdd : pos--; tab[pos-1] += tab[pos]; continue;
4226 case kSubstract : pos--; tab[pos-1] -= tab[pos]; continue;
4227 case kMultiply : pos--; tab[pos-1] *= tab[pos]; continue;
4228 case kDivide : pos--; if (tab[pos] == 0) tab[pos-1] = 0; // division by 0
4229 else tab[pos-1] /= tab[pos];
4230 continue;
4231 case kUnary : tab[pos] = (fPredefined[i]->fFunc10)(pdata[fOperOffset[i].fType0][fOperOffset[i].fOffset0]); pos++;continue;
4232 case kBinary : tab[pos] = (fPredefined[i]->fFunc110)(pdata[fOperOffset[i].fType0][fOperOffset[i].fOffset0],
4233 pdata[fOperOffset[i].fType1][fOperOffset[i].fOffset1]);pos++;continue;
4234
4235 case kThree : tab[pos] = (fPredefined[i]->fFunc1110)(pdata[fOperOffset[i].fType0][fOperOffset[i].fOffset0],
4236 pdata[fOperOffset[i].fType1][fOperOffset[i].fOffset1],
4237 pdata[fOperOffset[i].fType2][fOperOffset[i].fOffset2]); pos++; continue;
4238
4239 case kFDM : tab[pos] = (fPredefined[i]->fFuncG)(&x[fOperOffset[i].fType0],&params[fOperOffset[i].fOffset0]); pos++;continue;
4240 case kFD1 : tab[pos-1] =(fPredefined[i]->fFunc10)(tab[pos-1]); continue;
4241 case kFD2 : pos--; tab[pos-1] = (fPredefined[i]->fFunc110)(tab[pos-1],tab[pos]); continue;
4242 case kFD3 : pos-=2; tab[pos-1] = (fPredefined[i]->fFunc1110)(tab[pos-1],tab[pos],tab[pos+1]); continue;
4243 }
4244 //
4245 switch(opcode) {
4246 case kBoolOptimizeAnd: {
4247 if (!tab[pos-1]) i=fOperOffset[i].fToJump;
4248 continue;
4249 }
4250 case kBoolOptimizeOr: {
4251 if (tab[pos-1]) i=fOperOffset[i].fToJump;
4252 continue;
4253 }
4254 case kAnd : pos--; tab[pos-1] = (bool)tab[pos]; continue; // use the fact that other were check before - see bool optimize
4255 case kOr : pos--; tab[pos-1] = (bool)tab[pos]; continue;
4256 }
4257 switch(opcode) {
4258 // case kabs : tab[pos-1] = TMath::Abs(tab[pos-1]); continue;
4259 case kabs : if (tab[pos-1]<0) tab[pos-1]=-tab[pos-1]; continue;
4260 case ksign : if (tab[pos-1] < 0) tab[pos-1] = -1; else tab[pos-1] = 1; continue;
4261 case kint : tab[pos-1] = Double_t(Int_t(tab[pos-1])); continue;
4262 case kpow : pos--; tab[pos-1] = TMath::Power(tab[pos-1],tab[pos]); continue;
4263
4264 case kModulo : {pos--;
4265 Long64_t int1((Long64_t)tab[pos-1]);
4266 Long64_t int2((Long64_t)tab[pos]);
4267 tab[pos-1] = Double_t(int1%int2);
4268 continue;}
4269
4270
4271 case kStringConst: { strpos++; stringStack[strpos-1] = (char*)fExprOptimized[i].Data(); pos++; tab[pos-1] = 0; continue; }
4272 case kfmod : pos--; tab[pos-1] = fmod(tab[pos-1],tab[pos]); continue;
4273
4274 case kstrstr : strpos -= 2; pos-=2; pos++;
4275 if (strstr(stringStack[strpos],stringStack[strpos+1])) tab[pos-1]=1;
4276 else tab[pos-1]=0;
4277 continue;
4278 case kpi : pos++; tab[pos-1] = TMath::ACos(-1); continue;
4279
4280
4281 case kSignInv: tab[pos-1] = -1 * tab[pos-1]; continue;
4282
4283 case krndm : pos++; tab[pos-1] = gRandom->Rndm(); continue;
4284
4285
4286 case kEqual: pos--; if (tab[pos-1] == tab[pos]) tab[pos-1]=1;
4287 else tab[pos-1]=0;
4288 continue;
4289 case kNotEqual : pos--; if (tab[pos-1] != tab[pos]) tab[pos-1]=1;
4290 else tab[pos-1]=0;
4291 continue;
4292 case kNot : if (tab[pos-1]!=0) tab[pos-1] = 0; else tab[pos-1] = 1;
4293 continue;
4294
4295 case kStringEqual : strpos -= 2; pos -=2 ; pos++;
4296 if (!strcmp(stringStack[strpos+1],stringStack[strpos])) tab[pos-1]=1;
4297 else tab[pos-1]=0;
4298 continue;
4299 case kStringNotEqual: strpos -= 2; pos -= 2; pos++;
4300 if (strcmp(stringStack[strpos+1],stringStack[strpos])) tab[pos-1]=1;
4301 else tab[pos-1]=0;
4302 continue;
4303
4304 case kBitAnd : pos--; tab[pos-1]= ((Int_t) tab[pos-1]) & ((Int_t) tab[pos]); continue;
4305 case kBitOr : pos--; tab[pos-1]= ((Int_t) tab[pos-1]) | ((Int_t) tab[pos]); continue;
4306 case kLeftShift : pos--; tab[pos-1]= ((Int_t) tab[pos-1]) <<((Int_t) tab[pos]); continue;
4307 case kRightShift: pos--; tab[pos-1]= ((Int_t) tab[pos-1]) >>((Int_t) tab[pos]); continue;
4308
4309 case kJump : i = (oper & kTFOperMask); continue;
4310 case kJumpIf : pos--; if (!tab[pos]) i = (oper & kTFOperMask); continue;
4311
4312 case kBoolOptimize: {
4313 // boolean operation optimizer
4314
4315 int param = (oper & kTFOperMask);
4316 int op = param % 10; // 1 is && , 2 is ||
4317
4318 if (op == 1 && (!tab[pos-1]) ) {
4319 // &&: skip the right part if the left part is already false
4320
4321 i += param / 10;
4322
4323 // Preserve the existing behavior (i.e. the result of a&&b is
4324 // either 0 or 1)
4325 tab[pos-1] = 0;
4326
4327 } else if (op == 2 && tab[pos-1] ) {
4328 // ||: skip the right part if the left part is already true
4329
4330 i += param / 10;
4331
4332 // Preserve the existing behavior (i.e. the result of a||b is
4333 // either 0 or 1)
4334 tab[pos-1] = 1;
4335
4336 }
4337
4338 continue;
4339 }
4340
4341 }
4342 switch(opcode) {
4343
4344#define R__EXPO(var) \
4345 { \
4346 pos++; int param = (oper & kTFOperMask); \
4347 tab[pos-1] = TMath::Exp(params[param]+params[param+1]*x[var]); \
4348 continue; \
4349 }
4350 // case kexpo:
4351 case kxexpo: R__EXPO(0);
4352 case kyexpo: R__EXPO(1);
4353 case kzexpo: R__EXPO(2);
4354 case kxyexpo:{ pos++; int param = (oper & kTFOperMask);
4355 tab[pos-1] = TMath::Exp(params[param]+params[param+1]*x[0]+params[param+2]*x[1]);
4356 continue; }
4357#ifdef R__GAUS
4358#undef R__GAUS
4359#endif
4360#define R__GAUS(var) \
4361 { \
4362 pos++; int param = (oper & kTFOperMask); \
4363 tab[pos-1] = params[param]*TMath::Gaus(x[var],params[param+1], \
4364 params[param+2],IsNormalized()); \
4365 continue; \
4366 }
4367
4368 // case kgaus:
4369 case kxgaus: R__GAUS(0);
4370 case kygaus: R__GAUS(1);
4371 case kzgaus: R__GAUS(2);
4372 case kxygaus: { pos++; int param = (oper & kTFOperMask);
4373 Double_t intermede1;
4374 if (params[param+2] == 0) {
4375 intermede1=1e10;
4376 } else {
4377 intermede1=Double_t((x[0]-params[param+1])/params[param+2]);
4378 }
4379 Double_t intermede2;
4380 if (params[param+4] == 0) {
4381 intermede2=1e10;
4382 } else {
4383 intermede2=Double_t((x[1]-params[param+3])/params[param+4]);
4384 }
4385 tab[pos-1] = params[param]*TMath::Exp(-0.5*(intermede1*intermede1+intermede2*intermede2));
4386 continue; }
4387
4388#define R__LANDAU(var) \
4389 { \
4390 pos++; const int param = (oper & kTFOperMask); \
4391 tab[pos-1] = params[param]*TMath::Landau(x[var],params[param+1],params[param+2],IsNormalized()); \
4392 continue; \
4393 }
4394 // case klandau:
4395 case kxlandau: R__LANDAU(0);
4396 case kylandau: R__LANDAU(1);
4397 case kzlandau: R__LANDAU(2);
4398 case kxylandau: { pos++; int param = oper&0x7fffff /* ActionParams[i] */ ;
4399 Double_t intermede1=TMath::Landau(x[0], params[param+1], params[param+2],IsNormalized());
4400 Double_t intermede2=TMath::Landau(x[1], params[param+3], params[param+4],IsNormalized());
4401 tab[pos-1] = params[param]*intermede1*intermede2;
4402 continue;
4403 }
4404
4405#define R__POLY(var) \
4406 { \
4407 pos++; int param = (oper & kTFOperMask); \
4408 tab[pos-1] = 0; Double_t intermede = 1; \
4409 Int_t inter = param/100; /* arrondit */ \
4410 Int_t int1= param-inter*100-1; /* aucune simplification ! (sic) */ \
4411 for (j=0 ;j<inter+1;j++) { \
4412 tab[pos-1] += intermede*params[j+int1]; \
4413 intermede *= x[var]; \
4414 } \
4415 continue; \
4416 }
4417 // case kpol:
4418 case kxpol: R__POLY(0);
4419 case kypol: R__POLY(1);
4420 case kzpol: R__POLY(2);
4421
4422 case kDefinedVariable : {
4423 if (!precalculated) {
4424 precalculated = 1;
4425 for(j=0;j<fNval;j++) param_calc[j]=DefinedValue(j);
4426 }
4427 pos++; tab[pos-1] = param_calc[(oper & kTFOperMask)];
4428 continue;
4429 }
4430
4431 case kDefinedString : {
4432 int param = (oper & kTFOperMask);
4433 if (!precalculated_str) {
4434 precalculated_str=1;
4435 for (j=0;j<fNstring;j++) string_calc[j]=DefinedString(j);
4436 }
4437 strpos++; stringStack[strpos-1] = string_calc[param];
4438 pos++; tab[pos-1] = 0;
4439 continue;
4440 }
4441
4442 case kFunctionCall: {
4443 // an external function call
4444
4445 int param = (oper & kTFOperMask);
4446 int fno = param / 1000;
4447 int nargs = param % 1000;
4448
4449 // Retrieve the function
4450 TMethodCall *method = (TMethodCall*)fFunctions.At(fno);
4451
4452 // Set the arguments
4453 method->ResetParam();
4454 if (nargs) {
4455 UInt_t argloc = pos-nargs;
4456 for(j=0;j<nargs;j++,argloc++,pos--) {
4457 method->SetParam(tab[argloc]);
4458 }
4459 }
4460 pos++;
4461 Double_t ret;
4462 method->Execute(ret);
4463 tab[pos-1] = ret; // check for the correct conversion!
4464
4465 continue;
4466 };
4467 }
4468 if (!TestBit(kOptimizationError)) {
4469 SetBit(kOptimizationError);
4470 Warning("EvalParFast","Found an unsupported optmized opcode (%d)",oper >> kTFOperShift);
4471 }
4472 }
4473 Double_t result0 = tab[0];
4474 return result0;
4475
4476}
4477
4478////////////////////////////////////////////////////////////////////////////////
4479/// Pre compile function
4480
4481Int_t TFormula::PreCompile()
4482{
4483 TString str = fTitle;
4484 if (str.Length()<3) return 1;
4485 if (str[str.Length()-1]!='+'&&str[str.Length()-2]!='+') return 1;
4486 str[str.Length()-2]=0;
4487 TString funName("preformula_");
4488 funName += fName;
4489 if (ROOT::v5::TFormulaPrimitive::FindFormula(funName)) return 0;
4490 TString fileName;
4491 fileName.Form("/tmp/%s.C",funName.Data());
4492
4493 FILE *hf;
4494 hf = fopen(fileName.Data(),"w");
4495 if (hf == 0) {
4496 Error("PreCompile","Unable to open the file %s for writing.",fileName.Data());
4497 return 1;
4498 }
4499 fprintf(hf, "/////////////////////////////////////////////////////////////////////////\n");
4500 fprintf(hf, "// This code has been automatically generated \n");
4501 //
4502 fprintf(hf, "Double_t %s(Double_t *x, Double_t *p){",funName.Data());
4503 fprintf(hf, "return (%s);\n}",str.Data());
4504
4505 // fprintf("ROOT::v5::TFormulaPrimitive::AddFormula(new ROOT::v5::TFormulaPrimitive(\"%s::%s\",\"%s::%s\",(ROOT::v5::TFormulaPrimitive::GenFunc0)%s::%s));\n",
4506 // clname,method->GetName(),clname,method->GetName(),clname,method->GetName());
4507 fclose(hf);
4508
4509 return 0;
4510
4511
4512}
4513
4514////////////////////////////////////////////////////////////////////////////////
4515/// static function to set the maximum value of 3 parameters
4516///
4517/// - maxop : maximum number of operations
4518/// - maxpar : maximum number of parameters
4519/// - maxconst : maximum number of constants
4520///
4521/// None of these parameters cannot be less than 10 (default is 1000)
4522/// call this function to increase one or all maxima when processing
4523/// very complex formula, eg TFormula::SetMaxima(100000,1000,1000000);
4524/// If you process many functions with a small number of operations/parameters
4525/// you may gain some memory and performance by decreasing these values.
4526
4527void TFormula::SetMaxima(Int_t maxop, Int_t maxpar, Int_t maxconst)
4528{
4529 gMAXOP = TMath::Max(10,maxop);
4530 gMAXPAR = TMath::Max(10,maxpar);
4531 gMAXCONST = TMath::Max(10,maxconst);
4532}
4533
4534////////////////////////////////////////////////////////////////////////////////
4535/// static function to get the maximum value of 3 parameters
4536/// -maxop : maximum number of operations
4537/// -maxpar : maximum number of parameters
4538/// -maxconst : maximum number of constants
4539
4540void TFormula::GetMaxima(Int_t& maxop, Int_t& maxpar, Int_t& maxconst)
4541{
4542 maxop = gMAXOP;
4543 maxpar = gMAXPAR;
4544 maxconst = gMAXCONST;
4545}
4546
4547 } // end namespace v5
4548
4549} // end namespace ROOT
b
#define b(i)
Definition RSha256.hxx:100
f
#define f(i)
Definition RSha256.hxx:104
c
#define c(i)
Definition RSha256.hxx:101
s1
#define s1(x)
Definition RSha256.hxx:91
e
#define e(i)
Definition RSha256.hxx:103
Int_t
int Int_t
Definition RtypesCore.h:45
Version_t
short Version_t
Definition RtypesCore.h:65
UInt_t
unsigned int UInt_t
Definition RtypesCore.h:46
kFALSE
const Bool_t kFALSE
Definition RtypesCore.h:92
ULong_t
unsigned long ULong_t
Definition RtypesCore.h:55
Long_t
long Long_t
Definition RtypesCore.h:54
Double_t
double Double_t
Definition RtypesCore.h:59
Long64_t
long long Long64_t
Definition RtypesCore.h:73
Float_t
float Float_t
Definition RtypesCore.h:57
kTRUE
const Bool_t kTRUE
Definition RtypesCore.h:91
Option_t
const char Option_t
Definition RtypesCore.h:66
BIT
#define BIT(n)
Definition Rtypes.h:85
ClassImp
#define ClassImp(name)
Definition Rtypes.h:364
R__ASSERT
#define R__ASSERT(e)
Definition TError.h:120
R__GAUS
#define R__GAUS(var)
R__LANDAU
#define R__LANDAU(var)
gMAXCONST
static Int_t gMAXCONST
Definition TFormula_v5.cxx:34
gMAXSTRINGFOUND
const Int_t gMAXSTRINGFOUND
Definition TFormula_v5.cxx:35
R__POLY
#define R__POLY(var)
gMAXPAR
static Int_t gMAXPAR
Definition TFormula_v5.cxx:34
R__EXPO
#define R__EXPO(var)
gMAXOP
static Int_t gMAXOP
Definition TFormula_v5.cxx:34
kOptimizationError
const UInt_t kOptimizationError
Definition TFormula_v5.cxx:36
name
char name[80]
Definition TGX11.cxx:110
gInterpreter
#define gInterpreter
Definition TInterpreter.h:562
gROOTMutex
R__EXTERN TVirtualMutex * gROOTMutex
Definition TROOT.h:63
gROOT
#define gROOT
Definition TROOT.h:406
gRandom
R__EXTERN TRandom * gRandom
Definition TRandom.h:62
Form
char * Form(const char *fmt,...)
Printf
void Printf(const char *fmt,...)
R__LOCKGUARD
#define R__LOCKGUARD(mutex)
Definition TVirtualMutex.h:104
proto
const char * proto
Definition civetweb.c:16604
snprintf
#define snprintf
Definition civetweb.c:1540
ROOT::v5::TFormulaPrimitive
The Formula Primitive class.
Definition TFormulaPrimitive.h:38
ROOT::v5::TFormulaPrimitive::fFunc1110
GenFunc1110 fFunc1110
pointer to the function
Definition TFormulaPrimitive.h:59
ROOT::v5::TFormulaPrimitive::fFunc10
GenFunc10 fFunc10
pointer to the function
Definition TFormulaPrimitive.h:57
ROOT::v5::TFormulaPrimitive::FindFormula
static TFormulaPrimitive * FindFormula(const char *name)
Find the formula in the list of formulas.
Definition TFormulaPrimitive_v5.cxx:309
ROOT::v5::TFormulaPrimitive::TFuncG
Double_t(TObject::* TFuncG)(const Double_t *, const Double_t *) const
Definition TFormulaPrimitive.h:46
ROOT::v5::TFormulaPrimitive::fFunc110
GenFunc110 fFunc110
pointer to the function
Definition TFormulaPrimitive.h:58
ROOT::v5::TFormulaPrimitive::Eval
Double_t Eval(Double_t *x)
Eval primitive function at point x.
Definition TFormulaPrimitive_v5.cxx:209
ROOT::v5::TFormula
The FORMULA class (ROOT version 5)
Definition TFormula.h:65
ROOT::v5::TFormula::fConst
Double_t * fConst
Definition TFormula.h:82
ROOT::v5::TFormula::GetNumber
virtual Int_t GetNumber() const
Definition TFormula.h:239
ROOT::v5::TFormula::StringToNumber
virtual Bool_t StringToNumber(Int_t code)
Try to 'demote' a string into an array bytes.
Definition TFormula_v5.cxx:2202
ROOT::v5::TFormula::GetParameter
Double_t GetParameter(Int_t ipar) const
Return value of parameter number ipar.
Definition TFormula_v5.cxx:3222
ROOT::v5::TFormula::fExprOptimized
TString * fExprOptimized
Number of operators after optimization.
Definition TFormula.h:92
ROOT::v5::TFormula::~TFormula
virtual ~TFormula()
Formula default destructor.
Definition TFormula_v5.cxx:314
ROOT::v5::TFormula::IsNormalized
virtual Bool_t IsNormalized() const
Definition TFormula.h:248
ROOT::v5::TFormula::DefinedValue
virtual Double_t DefinedValue(Int_t code)
Return value corresponding to special code.
Definition TFormula_v5.cxx:2586
ROOT::v5::TFormula::Eval
virtual Double_t Eval(Double_t x, Double_t y=0, Double_t z=0, Double_t t=0) const
Evaluate this formula.
Definition TFormula_v5.cxx:2668
ROOT::v5::TFormula::fFunctions
TObjArray fFunctions
Definition TFormula.h:85
ROOT::v5::TFormula::GetParNumber
virtual Int_t GetParNumber(const char *name) const
Return parameter number by name.
Definition TFormula_v5.cxx:3255
ROOT::v5::TFormula::Convert
virtual void Convert(UInt_t fromVersion)
Definition TFormula_v5.cxx:3564
ROOT::v5::TFormula::SetActionOptimized
void SetActionOptimized(Int_t code, Int_t value, Int_t param=0)
Definition TFormula.h:115
ROOT::v5::TFormula::EvalPrimitive1
Double_t EvalPrimitive1(const Double_t *x, const Double_t *params)
Evaluate primitive formula.
Definition TFormula_v5.cxx:4144
ROOT::v5::TFormula::SetParName
virtual void SetParName(Int_t ipar, const char *name)
Set name of parameter number ipar.
Definition TFormula_v5.cxx:3432
ROOT::v5::TFormula::SetParameter
virtual void SetParameter(const char *name, Double_t parvalue)
Initialize parameter number ipar.
Definition TFormula_v5.cxx:3375
ROOT::v5::TFormula::EvalPrimitive0
Double_t EvalPrimitive0(const Double_t *x, const Double_t *params)
Evaluate primitive formula.
Definition TFormula_v5.cxx:4135
ROOT::v5::TFormula::Update
virtual void Update()
Definition TFormula.h:263
ROOT::v5::TFormula::AnalyzeFunction
virtual Bool_t AnalyzeFunction(TString &chaine, Int_t &err, Int_t offset=0)
Check if the chain as function call.
Definition TFormula_v5.cxx:347
ROOT::v5::TFormula::TFormula
TFormula()
Formula default constructor.
Definition TFormula_v5.cxx:136
ROOT::v5::TFormula::EvalPrimitive2
Double_t EvalPrimitive2(const Double_t *x, const Double_t *params)
Evaluate primitive formula.
Definition TFormula_v5.cxx:4153
ROOT::v5::TFormula::fExpr
TString * fExpr
Definition TFormula.h:78
ROOT::v5::TFormula::EvalParOld
virtual Double_t EvalParOld(const Double_t *x, const Double_t *params=0)
Evaluate this formula.
Definition TFormula_v5.cxx:2687
ROOT::v5::TFormula::fLinearParts
TObjArray fLinearParts
Definition TFormula.h:86
ROOT::v5::TFormula::ClearFormula
void ClearFormula(Option_t *option="")
Resets the objects.
Definition TFormula_v5.cxx:2223
ROOT::v5::TFormula::GetParName
virtual const char * GetParName(Int_t ipar) const
Return name of one parameter.
Definition TFormula_v5.cxx:3245
ROOT::v5::TFormula::SetNumber
virtual void SetNumber(Int_t number)
Definition TFormula.h:251
ROOT::v5::TFormula::fOperOffset
TOperOffset * fOperOffset
[fNOperOptimized] List of operators. (See documentation for changes made at version 7)
Definition TFormula.h:94
ROOT::v5::TFormula::DefinedString
virtual char * DefinedString(Int_t code)
Return address of string corresponding to special code.
Definition TFormula_v5.cxx:2555
ROOT::v5::TFormula::Optimize
void Optimize()
MI include.
De