ROOT  6.06/09
Reference Guide
TFormulaPrimitive_v5.cxx
Go to the documentation of this file.
1 // @(#)root/hist:$Id$
2 // Author: Marian Ivanov, 2005
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 <math.h>
13 
14 #include "v5/TFormulaPrimitive.h"
15 #include "TMath.h"
16 #include "TVirtualMutex.h"
17 #ifdef WIN32
18 #pragma optimize("",off)
19 #endif
20 
21 static TVirtualMutex* gTFormulaPrimativeListMutex = 0;
22 
23 
25 
26 namespace ROOT {
27  namespace v5 {
28 
29  void TMath_GenerInterface();
30 
31 
32 //______________________________________________________________________________
33 // The Formula Primitive class
34 //
35 // Helper class for TFormula to speed up TFormula evaluation
36 // TFormula can use all functions registered in the list of TFormulaPrimitives
37 // User can add new function to the list of primitives
38 // if FormulaPrimitive with given name is already defined new primitive is ignored
39 // Example:
40 // TFormulaPrimitive::AddFormula(new TFormulaPrimitive("Pow2","Pow2",TFastFun::Pow2));
41 // TF1 f1("f1","Pow2(x)");
42 //
43 //
44 //
45 // TFormulaPrimitive is used to get direct acces to the function pointers
46 // GenFunc - pointers to the static function
47 // TFunc - pointers to the data member functions
48 //
49 // The following sufixes are currently used, to describe function arguments:
50 // ------------------------------------------------------------------------
51 // G - generic layout - pointer to double (arguments), pointer to double (parameters)
52 // 10 - double
53 // 110 - double, double
54 // 1110 - double, double, double
55 
56 
57 //______________________________________________________________________________
58 // TFormula primitive
59 //
60 TObjArray * TFormulaPrimitive::fgListOfFunction = 0;
61 #ifdef R__COMPLETE_MEM_TERMINATION
62 namespace {
63  class TFormulaPrimitiveCleanup {
64  TObjArray **fListOfFunctions;
65  public:
66  TFormulaPrimitiveCleanup(TObjArray **functions) : fListOfFunctions(functions) {}
67  ~TFormulaPrimitiveCleanup() {
68  delete *fListOfFunctions;
69  }
70  };
71 
72 }
73 #endif
74 
75 ////////////////////////////////////////////////////////////////////////////////
76 /// Default constructor.
77 
78 TFormulaPrimitive::TFormulaPrimitive() : TNamed(),
79  fFuncG(0),
80  fType(0),fNArguments(0),fNParameters(0),fIsStatic(kTRUE)
81 {
82 }
83 
84 
85 ////////////////////////////////////////////////////////////////////////////////
86 /// Constructor.
87 
88 TFormulaPrimitive::TFormulaPrimitive(const char *name,const char *formula,
89  GenFunc0 fpointer) : TNamed(name,formula),
90  fFunc0(fpointer),
91  fType(0),fNArguments(0),fNParameters(0),fIsStatic(kTRUE)
92 {
93 }
94 
95 
96 ////////////////////////////////////////////////////////////////////////////////
97 /// Constructor.
98 
99 TFormulaPrimitive::TFormulaPrimitive(const char *name,const char *formula,
100  GenFunc10 fpointer) : TNamed(name,formula),
101  fFunc10(fpointer),
102  fType(10),fNArguments(1),fNParameters(0),fIsStatic(kTRUE)
103 {
104 }
105 
106 
107 ////////////////////////////////////////////////////////////////////////////////
108 /// Constructor.
109 
110 TFormulaPrimitive::TFormulaPrimitive(const char *name,const char *formula,
111  GenFunc110 fpointer) : TNamed(name,formula),
112  fFunc110(fpointer),
113  fType(110),fNArguments(2),fNParameters(0),fIsStatic(kTRUE)
114 {
115 }
116 
117 
118 ////////////////////////////////////////////////////////////////////////////////
119 /// Constructor.
120 
121 TFormulaPrimitive::TFormulaPrimitive(const char *name,const char *formula,
122  GenFunc1110 fpointer) : TNamed(name,formula),
123  fFunc1110(fpointer),
124  fType(1110),fNArguments(3),fNParameters(0),fIsStatic(kTRUE)
125 {
126 }
127 
128 
129 ////////////////////////////////////////////////////////////////////////////////
130 /// Constructor.
131 
132 TFormulaPrimitive::TFormulaPrimitive(const char *name,const char *formula,
133  GenFuncG fpointer,Int_t npar) : TNamed(name,formula),
134  fFuncG(fpointer),
135  fType(-1),fNArguments(2),fNParameters(npar),fIsStatic(kTRUE)
136 {
137 }
138 
139 
140 ////////////////////////////////////////////////////////////////////////////////
141 /// Constructor.
142 
143 TFormulaPrimitive::TFormulaPrimitive(const char *name,const char *formula,
144  TFuncG fpointer) : TNamed(name,formula),
145  fTFuncG(fpointer),
146  fType(0),fNArguments(0),fNParameters(0),fIsStatic(kFALSE)
147 {
148 }
149 
150 
151 ////////////////////////////////////////////////////////////////////////////////
152 /// Constructor.
153 
154 TFormulaPrimitive::TFormulaPrimitive(const char *name,const char *formula,
155  TFunc0 fpointer) : TNamed(name,formula),
156  fTFunc0(fpointer),
157  fType(0),fNArguments(0),fNParameters(0),fIsStatic(kFALSE)
158 {
159 }
160 
161 
162 ////////////////////////////////////////////////////////////////////////////////
163 /// Constructor.
164 
165 TFormulaPrimitive::TFormulaPrimitive(const char *name,const char *formula,
166  TFunc10 fpointer) : TNamed(name,formula),
167  fTFunc10(fpointer),
168  fType(-10),fNArguments(1),fNParameters(0),fIsStatic(kFALSE)
169 {
170 }
171 
172 
173 ////////////////////////////////////////////////////////////////////////////////
174 /// Constructor.
175 
176 TFormulaPrimitive::TFormulaPrimitive(const char *name,const char *formula,
177  TFunc110 fpointer) : TNamed(name,formula),
178  fTFunc110(fpointer),
179  fType(-110),fNArguments(2),fNParameters(0),fIsStatic(kFALSE)
180 {
181 }
182 
183 
184 ////////////////////////////////////////////////////////////////////////////////
185 /// Constructor.
186 
187 TFormulaPrimitive::TFormulaPrimitive(const char *name,const char *formula,
188  TFunc1110 fpointer) :TNamed(name,formula),
189  fTFunc1110(fpointer),
190  fType(-1110),fNArguments(3),fNParameters(0),fIsStatic(kFALSE)
191 {
192 }
193 
194 
195 ////////////////////////////////////////////////////////////////////////////////
196 /// Add formula to the list of primitive formulas.
197 /// If primitive formula already defined do nothing.
198 
199 Int_t TFormulaPrimitive::AddFormula(TFormulaPrimitive * formula)
200 {
201  R__LOCKGUARD2(gTFormulaPrimativeListMutex);
202  if (fgListOfFunction == 0) BuildBasicFormulas();
203  if (FindFormula(formula->GetName(),formula->fNArguments)){
204  delete formula;
205  return 0;
206  }
207  fgListOfFunction->AddLast(formula);
208  return 1;
209 }
210 
211 
212 ////////////////////////////////////////////////////////////////////////////////
213 /// Eval primitive function at point x.
214 
216 {
217  if (fIsStatic == kFALSE) return 0;
218 
219  if (fType==0) return fFunc0();
220  if (fType==10) {
221  return fFunc10(x[0]);
222  }
223  if (fType==110) {
224  return fFunc110(x[0],x[1]);
225  }
226  if (fType==1110) {
227  return fFunc1110(x[0],x[1],x[2]);
228  }
229  return 0;
230 }
231 
232 
233 ////////////////////////////////////////////////////////////////////////////////
234 /// Eval member function of object o at point x.
235 
237 {
238  if (fIsStatic == kTRUE) return 0;
239  if (fType== 0) return (*o.*fTFunc0)();
240  if (fType==-10) return (*o.*fTFunc10)(*x);
241  if (fType==-110) return (*o.*fTFunc110)(x[0],x[1]);
242  if (fType==-1110) return (*o.*fTFunc1110)(x[0],x[1],x[2]);
243  return 0;
244 }
245 
246 
247 ////////////////////////////////////////////////////////////////////////////////
248 /// Eval primitive parametric function.
249 
251 {
252  return fFuncG(x,param);
253 }
254 
255 
256 #define RTFastFun__POLY(var) \
257 { \
258  Double_t res= param[var-1]+param[var]*x[0]; \
259  for (Int_t j=var-1 ;j>0;j--) res = param[j-1]+x[0]*res; \
260  return res; \
261 }
262 
263 namespace TFastFun {
264  //
265  // Namespace with basic primitive functions registered by TFormulaPrimitive
266  // all function registered by TFormulaPrimitive can be used in TFormula
267  //
268  Double_t Pow2(Double_t x){return x*x;}
269  Double_t Pow3(Double_t x){return x*x*x;}
270  Double_t Pow4(Double_t x){return x*x*x*x;}
271  Double_t Pow5(Double_t x){return x*x*x*x*x;}
272  inline Double_t FPoln(Double_t *x, Double_t *param, Int_t npar);
273  Double_t FPol0(Double_t * /*x*/, Double_t *param){ return param[0];}
274  Double_t FPol1(Double_t *x, Double_t *param){ return param[0]+param[1]*x[0];}
275  Double_t FPol2(Double_t *x, Double_t *param){ return param[0]+x[0]*(param[1]+param[2]*x[0]);}
276  Double_t FPol3(Double_t *x, Double_t *param){ return param[0]+x[0]*(param[1]+x[0]*(param[2]+param[3]*x[0]));}
277  Double_t FPol4(Double_t *x, Double_t *param){ RTFastFun__POLY(4)}
278  Double_t FPol5(Double_t *x, Double_t *param){ RTFastFun__POLY(5)}
279  Double_t FPol6(Double_t *x, Double_t *param){ RTFastFun__POLY(6)}
280  Double_t FPol7(Double_t *x, Double_t *param){ RTFastFun__POLY(7)}
281  Double_t FPol8(Double_t *x, Double_t *param){ RTFastFun__POLY(8)}
282  Double_t FPol9(Double_t *x, Double_t *param){ RTFastFun__POLY(9)}
283  Double_t FPol10(Double_t *x, Double_t *param){ RTFastFun__POLY(10)}
284  //
285  //
286  Double_t PlusXY(Double_t x,Double_t y){return x+y;}
287  Double_t MinusXY(Double_t x,Double_t y){return x-y;}
288  Double_t MultXY(Double_t x,Double_t y){return x*y;}
289  Double_t DivXY(Double_t x, Double_t y){return TMath::Abs(y)>0 ? x/y:0;}
290  Double_t XpYpZ(Double_t x, Double_t y, Double_t z){ return x+y+z;}
291  Double_t XxYxZ(Double_t x, Double_t y, Double_t z){ return x*y*z;}
292  Double_t XxYpZ(Double_t x, Double_t y, Double_t z){ return x*(y+z);}
293  Double_t XpYxZ(Double_t x, Double_t y, Double_t z){ return x+(y*z);}
294  Double_t Gaus(Double_t x, Double_t mean, Double_t sigma);
295  Double_t Gausn(Double_t x, Double_t mean, Double_t sigma);
296  Double_t Landau(Double_t x, Double_t mean, Double_t sigma){return TMath::Landau(x,mean,sigma,kFALSE);}
297  Double_t Landaun(Double_t x, Double_t mean, Double_t sigma){return TMath::Landau(x,mean,sigma,kTRUE);}
298  Double_t Sqrt(Double_t x) {return x>0?sqrt(x):0;}
299  //
300  Double_t Sign(Double_t x){return (x<0)? -1:1;}
301  Double_t Nint(Double_t x){return TMath::Nint(x);}
302  Double_t Abs(Double_t x){return TMath::Abs(x);}
303  //logical
304  Double_t XandY(Double_t x, Double_t y){ return (x*y>0.1);}
305  Double_t XorY(Double_t x, Double_t y) { return (x+y>0.1);}
306  Double_t XgY(Double_t x, Double_t y) {return (x>y);}
307  Double_t XgeY(Double_t x, Double_t y) {return (x>=y);}
308  Double_t XlY(Double_t x, Double_t y) {return (x<y);}
309  Double_t XleY(Double_t x, Double_t y) {return (x<=y);}
310  Double_t XeY(Double_t x,Double_t y) {return (x==y);}
311  Double_t XneY(Double_t x,Double_t y) {return (x!=y);}
312  Double_t XNot(Double_t x){ return (x<0.1);}
313 };
314 
315 
316 ////////////////////////////////////////////////////////////////////////////////
317 /// Find the formula in the list of formulas.
318 
319 TFormulaPrimitive* TFormulaPrimitive::FindFormula(const char* name)
320 {
321  R__LOCKGUARD2(gTFormulaPrimativeListMutex);
322  if (!fgListOfFunction) {
323  BuildBasicFormulas();
324  }
325  Int_t nobjects = fgListOfFunction->GetEntries();
326  for (Int_t i = 0; i < nobjects; ++i) {
327  TFormulaPrimitive *formula = (TFormulaPrimitive*)fgListOfFunction->At(i);
328  if (formula && 0==strcmp(name, formula->GetName())) return formula;
329  }
330  return 0;
331 }
332 
333 
334 ////////////////////////////////////////////////////////////////////////////////
335 /// Find the formula in the list of formulas.
336 
337 TFormulaPrimitive* TFormulaPrimitive::FindFormula(const char* name, UInt_t nargs)
338 {
339  R__LOCKGUARD2(gTFormulaPrimativeListMutex);
340  if (!fgListOfFunction) {
341  BuildBasicFormulas();
342  }
343  Int_t nobjects = fgListOfFunction->GetEntries();
344  for (Int_t i = 0; i < nobjects; ++i) {
345  TFormulaPrimitive *prim = (TFormulaPrimitive*)fgListOfFunction->At(i);
346  if (prim) {
347  bool match = ( ((UInt_t)prim->fNArguments) == nargs );
348  if (match && 0==strcmp(name, prim->GetName())) return prim;
349  }
350  }
351  return 0;
352 }
353 
354 
355 ////////////////////////////////////////////////////////////////////////////////
356 /// Find the formula in the list of formulas.
357 
358 TFormulaPrimitive* TFormulaPrimitive::FindFormula(const char* name, const char *args)
359 {
360  // let's count the argument(s)
361  if (args) {
362  Int_t nargs = 0;
363  if (args[0]!=')') {
364  nargs = 1;
365  int nest = 0;
366  for(UInt_t c = 0; c < strlen(args); ++c ) {
367  switch (args[c]) {
368  case '(': ++nest; break;
369  case ')': --nest; break;
370  case '<': ++nest; break;
371  case '>': --nest; break;
372  case ',': nargs += (nest==0); break;
373  }
374  }
375  }
376  return FindFormula(name,nargs);
377  } else {
378  return FindFormula(name);
379  }
380  return 0;
381 }
382 
383 
384 ////////////////////////////////////////////////////////////////////////////////
385 /// FPoln.
386 
387 Double_t TFastFun::FPoln(Double_t *x, Double_t *param, Int_t npar)
388 {
389  Double_t res = 0; Double_t temp=1;
390  for (Int_t j=npar ;j>=0;j--) {
391  res += temp*param[j];
392  temp *= *x;
393  }
394  return res;
395 }
396 
397 
398 ////////////////////////////////////////////////////////////////////////////////
399 /// Gauss.
400 
402 {
403  if (sigma == 0) return 1.e30;
404  Double_t arg = (x-mean)/sigma;
405  return TMath::Exp(-0.5*arg*arg);
406 }
407 
408 
409 ////////////////////////////////////////////////////////////////////////////////
410 /// Normalize gauss.
411 
412 Double_t TFastFun::Gausn(Double_t x, Double_t mean, Double_t sigma)
413 {
414  if (sigma == 0) return 0;
415  Double_t arg = (x-mean)/sigma;
416  return TMath::Exp(-0.5*arg*arg)/(2.50662827463100024*sigma); //sqrt(2*Pi)=2.50662827463100024
417 }
418 
419 
420 ////////////////////////////////////////////////////////////////////////////////
421 /// Built-in functions.
422 
423 Int_t TFormulaPrimitive::BuildBasicFormulas()
424 {
425  R__LOCKGUARD2(gTFormulaPrimativeListMutex);
426  if (fgListOfFunction==0) {
427  fgListOfFunction = new TObjArray(1000);
428  fgListOfFunction->SetOwner(kTRUE);
429  }
430 #ifdef R__COMPLETE_MEM_TERMINATION
431  static TFormulaPrimitiveCleanup gCleanup(&fgListOfFunction);
432 #endif
433 
434  //
435  // logical
436  //
437  AddFormula(new TFormulaPrimitive("XandY","XandY",TFastFun::XandY));
438  AddFormula(new TFormulaPrimitive("XorY","XorY",TFastFun::XorY));
439  AddFormula(new TFormulaPrimitive("XNot","XNot",TFastFun::XNot));
440  AddFormula(new TFormulaPrimitive("XlY","XlY",TFastFun::XlY));
441  AddFormula(new TFormulaPrimitive("XleY","XleY",TFastFun::XleY));
442  AddFormula(new TFormulaPrimitive("XgY","XgY",TFastFun::XgY));
443  AddFormula(new TFormulaPrimitive("XgeY","XgeY",TFastFun::XgeY));
444  AddFormula(new TFormulaPrimitive("XeY","XeY",TFastFun::XeY));
445  AddFormula(new TFormulaPrimitive("XneY","XneY",TFastFun::XneY));
446  // addition + multiplication
447  AddFormula(new TFormulaPrimitive("PlusXY","PlusXY",TFastFun::PlusXY));
448  AddFormula(new TFormulaPrimitive("MinusXY","MinusXY",TFastFun::MinusXY));
449  AddFormula(new TFormulaPrimitive("MultXY","MultXY",TFastFun::MultXY));
450  AddFormula(new TFormulaPrimitive("DivXY","DivXY",TFastFun::DivXY));
451  AddFormula(new TFormulaPrimitive("XpYpZ","XpYpZ",TFastFun::XpYpZ));
452  AddFormula(new TFormulaPrimitive("XxYxZ","XxYxZ",TFastFun::XxYxZ));
453  AddFormula(new TFormulaPrimitive("XxYpZ","XxYpZ",TFastFun::XxYpZ));
454  AddFormula(new TFormulaPrimitive("XpYxZ","XpYxZ",TFastFun::XpYxZ));
455  //
456  //
457  AddFormula(new TFormulaPrimitive("Gaus","Gaus",TFastFun::Gaus));
458  AddFormula(new TFormulaPrimitive("Gausn","Gausn",TFastFun::Gausn));
459  AddFormula(new TFormulaPrimitive("Landau","Landau",TFastFun::Landau));
460  AddFormula(new TFormulaPrimitive("Landaun","Landaun",TFastFun::Landaun));
461  //
462  //
463  // polynoms
464  //
465  //
466  AddFormula(new TFormulaPrimitive("Pol0","Pol0",(GenFuncG)TFastFun::FPol0,1));
467  AddFormula(new TFormulaPrimitive("Pol1","Pol1",(GenFuncG)TFastFun::FPol1,2));
468  AddFormula(new TFormulaPrimitive("Pol2","Pol2",(GenFuncG)TFastFun::FPol2,3));
469  AddFormula(new TFormulaPrimitive("Pol3","Pol3",(GenFuncG)TFastFun::FPol3,4));
470  AddFormula(new TFormulaPrimitive("Pol4","Pol4",(GenFuncG)TFastFun::FPol4,5));
471  AddFormula(new TFormulaPrimitive("Pol5","Pol5",(GenFuncG)TFastFun::FPol5,6));
472  AddFormula(new TFormulaPrimitive("Pol6","Pol6",(GenFuncG)TFastFun::FPol6,7));
473  AddFormula(new TFormulaPrimitive("Pol7","Pol7",(GenFuncG)TFastFun::FPol7,8));
474  AddFormula(new TFormulaPrimitive("Pol8","Pol8",(GenFuncG)TFastFun::FPol8,9));
475  AddFormula(new TFormulaPrimitive("Pol9","Pol9",(GenFuncG)TFastFun::FPol9,10));
476  AddFormula(new TFormulaPrimitive("Pol10","Pol10",(GenFuncG)TFastFun::FPol10,11));
477  //
478  // pows
479  AddFormula(new TFormulaPrimitive("Pow2","Pow2",TFastFun::Pow2));
480  AddFormula(new TFormulaPrimitive("Pow3","Pow3",TFastFun::Pow3));
481  AddFormula(new TFormulaPrimitive("Pow4","Pow4",TFastFun::Pow4));
482  AddFormula(new TFormulaPrimitive("Pow5","Pow5",TFastFun::Pow5));
483  //
484  //
485  AddFormula(new TFormulaPrimitive("TMath::Cos","TMath::Cos",cos)); // 10
486  AddFormula(new TFormulaPrimitive("cos","cos",cos)); // 10
487  AddFormula(new TFormulaPrimitive("TMath::Sin","TMath::Sin",sin)); // 11
488  AddFormula(new TFormulaPrimitive("sin","sin",sin)); // 11
489  AddFormula(new TFormulaPrimitive("TMath::Tan","TMath::Tan",tan)); // 12
490  AddFormula(new TFormulaPrimitive("tan","tan",tan)); // 12
491  AddFormula(new TFormulaPrimitive("TMath::ACos","TMath::ACos",acos)); // 13
492  AddFormula(new TFormulaPrimitive("acos","acos",acos)); // 13
493  AddFormula(new TFormulaPrimitive("TMath::ASin","TMath::ASin",asin)); // 14
494  AddFormula(new TFormulaPrimitive("asin","asin",asin)); // 14
495  AddFormula(new TFormulaPrimitive("TMath::ATan","TMath::ATan",atan)); // 15
496  AddFormula(new TFormulaPrimitive("atan","atan",atan)); // 15
497  AddFormula(new TFormulaPrimitive("TMath::ATan2","TMath::ATan2",atan2)); // 16
498  AddFormula(new TFormulaPrimitive("atan2","atan2",atan2)); // 16
499  // kpow = 20, ksq = 21, ksqrt = 22,
500  AddFormula(new TFormulaPrimitive("pow","pow",TMath::Power)); //20
501  AddFormula(new TFormulaPrimitive("sq","sq",TFastFun::Pow2)); //21
502  AddFormula(new TFormulaPrimitive("sqrt","sqrt",TFastFun::Sqrt)); //22
503  // kmin = 24, kmax = 25,
504  AddFormula(new TFormulaPrimitive("min","min",(GenFunc110)TMath::Min)); //24
505  AddFormula(new TFormulaPrimitive("max","max",(GenFunc110)TMath::Max)); //25
506  // klog = 30, kexp = 31, klog10 = 32,
507  AddFormula(new TFormulaPrimitive("log","log",TMath::Log)); //30
508  AddFormula(new TFormulaPrimitive("exp","exp",TMath::Exp)); //31
509  AddFormula(new TFormulaPrimitive("log10","log10",TMath::Log10)); //32
510  //
511  // cosh 70 acosh 73
512  // sinh 71 asinh 74
513  // tanh 72 atanh 75
514  //
515  AddFormula(new TFormulaPrimitive("TMath::CosH","TMath::Cosh",cosh)); // 70
516  AddFormula(new TFormulaPrimitive("cosh","cosh",cosh)); // 70
517  AddFormula(new TFormulaPrimitive("TMath::SinH","TMath::SinH",sinh)); // 71
518  AddFormula(new TFormulaPrimitive("sinh","sinh",sinh)); // 71
519  AddFormula(new TFormulaPrimitive("TMath::TanH","TMath::Tanh",tanh)); // 72
520  AddFormula(new TFormulaPrimitive("tanh","tanh",tanh)); // 72
521  AddFormula(new TFormulaPrimitive("TMath::ACosH","TMath::ACosh",TMath::ACosH)); // 73
522  AddFormula(new TFormulaPrimitive("acosh","acosH",TMath::ACosH)); // 73
523  AddFormula(new TFormulaPrimitive("TMath::ASinH","TMath::ASinh",TMath::ASinH)); // 74
524  AddFormula(new TFormulaPrimitive("acosh","acosH",TMath::ASinH)); // 74
525  AddFormula(new TFormulaPrimitive("TMath::ATanH","TMath::ATanh",TMath::ATanH)); // 75
526  AddFormula(new TFormulaPrimitive("atanh","atanh",TMath::ATanH)); // 75
527  //
528  AddFormula(new TFormulaPrimitive("TMath::Abs","TMath::Abs",TMath::Abs));
529  AddFormula(new TFormulaPrimitive("TMath::BreitWigner","TMath::BreitWigner",TMath::BreitWigner));
530 
531  //Disable direct access to TMath::Landau for now because of the default parameter.
532  //AddFormula(new TFormulaPrimitive("TMath::Landau","TMath::Landau",(TFormulaPrimitive::GenFunc1110)TMath::Landau));
533 
535  return 1;
536 }
537 
538  } // end namespace v5
539 
540 } // end namespace ROOT
541 
Double_t ACosH(Double_t)
Definition: TMath.cxx:80
Double_t Landau(Double_t x, Double_t mpv=0, Double_t sigma=1, Bool_t norm=kFALSE)
The LANDAU function.
Definition: TMath.cxx:473
An array of TObjects.
Definition: TObjArray.h:39
Double_t BreitWigner(Double_t x, Double_t mean=0, Double_t gamma=1)
Calculate a Breit Wigner function with mean and gamma.
Definition: TMath.cxx:442
double tanh(double)
T1 Sign(T1 a, T2 b)
Definition: TMathBase.h:155
Namespace for new ROOT classes and functions.
Definition: ROOT.py:1
Double_t Log(Double_t x)
Definition: TMath.h:526
This class implements a mutex interface.
Definition: TVirtualMutex.h:34
Short_t Min(Short_t a, Short_t b)
Definition: TMathBase.h:170
int Int_t
Definition: RtypesCore.h:41
const Bool_t kFALSE
Definition: Rtypes.h:92
Short_t Abs(Short_t d)
Definition: TMathBase.h:110
double cos(double)
LongDouble_t Power(LongDouble_t x, LongDouble_t y)
Definition: TMath.h:501
ClassImp(ROOT::v5::TFormulaPrimitive) namespace ROOT
double sqrt(double)
Double_t x[n]
Definition: legend1.C:17
double acos(double)
The TNamed class is the base class for all named ROOT classes.
Definition: TNamed.h:33
double sinh(double)
UChar_t mod R__LOCKGUARD2(gSrvAuthenticateMutex)
Double_t ASinH(Double_t)
Definition: TMath.cxx:67
Double_t Log10(Double_t x)
Definition: TMath.h:529
double sin(double)
double cosh(double)
unsigned int UInt_t
Definition: RtypesCore.h:42
double asin(double)
Double_t Gaus(Double_t x, Double_t mean=0, Double_t sigma=1, Bool_t norm=kFALSE)
Calculate a gaussian function with mean and sigma.
Definition: TMath.cxx:453
PyObject * fType
#define RTFastFun__POLY(var)
Double_t Exp(Double_t x)
Definition: TMath.h:495
double Double_t
Definition: RtypesCore.h:55
double atan2(double, double)
Double_t y[n]
Definition: legend1.C:17
double atan(double)
#define name(a, b)
Definition: linkTestLib0.cpp:5
Mother of all ROOT objects.
Definition: TObject.h:58
Short_t Max(Short_t a, Short_t b)
Definition: TMathBase.h:202
double tan(double)
Double_t Sqrt(Double_t x)
Definition: TMath.h:464
const Bool_t kTRUE
Definition: Rtypes.h:91
Int_t Nint(T x)
Definition: TMath.h:480
Double_t ATanH(Double_t)
Definition: TMath.cxx:93
void TMath_GenerInterface()