CLHEP VERSION Reference Documentation
   
CLHEP Home Page     CLHEP Documentation     CLHEP Bug Reports

Landau.cc
Go to the documentation of this file.
1 // -*- C++ -*-
2 // $Id: Landau.cc,v 1.4 2003/10/10 17:40:39 garren Exp $
3 // ---------------------------------------------------------------------------
4 
7 #include <cmath>
8 #include <assert.h>
9 
10 using namespace std;
11 
12 namespace Genfun {
13 FUNCTION_OBJECT_IMP(Landau)
14 
16  _peak("Peak", 5.0 ,0, 10),
17  _width("Width",1.0,0, 10)
18 {}
19 
20 Landau::~Landau() {
21 }
22 
23 Landau::Landau(const Landau & right):
24 AbsFunction(right),
25 _peak(right._peak),
26 _width(right._width)
27 {
28 }
29 
30 double Landau::operator() (double x) const {
31  double s = _width.getValue();
32  double x0 = _peak.getValue();
33  double xs = x0 + 0.222782*s;
34  return _denlan((x-xs)/s)/s;
35 }
36 
38  return _peak;
39 }
40 
42  return _width;
43 }
44 
45 const Parameter & Landau::peak() const {
46  return _peak;
47 }
48 
49 const Parameter & Landau::width() const {
50  return _width;
51 }
52 
53 double Landau::_denlan(double x) const {
54  /* Initialized data */
55 
56  static float p1[5] = { (float).4259894875,(float)-.124976255,(float)
57  .039842437,(float)-.006298287635,(float).001511162253 };
58  static float q5[5] = { (float)1.,(float)156.9424537,(float)3745.310488,(
59  float)9834.698876,(float)66924.28357 };
60  static float p6[5] = { (float)1.000827619,(float)664.9143136,(float)
61  62972.92665,(float)475554.6998,(float)-5743609.109 };
62  static float q6[5] = { (float)1.,(float)651.4101098,(float)56974.73333,(
63  float)165917.4725,(float)-2815759.939 };
64  static float a1[3] = { (float).04166666667,(float)-.01996527778,(float)
65  .02709538966 };
66  static float a2[2] = { (float)-1.84556867,(float)-4.284640743 };
67  static float q1[5] = { (float)1.,(float)-.3388260629,(float).09594393323,(
68  float)-.01608042283,(float).003778942063 };
69  static float p2[5] = { (float).1788541609,(float).1173957403,(float)
70  .01488850518,(float)-.001394989411,(float)1.283617211e-4 };
71  static float q2[5] = { (float)1.,(float).7428795082,(float).3153932961,(
72  float).06694219548,(float).008790609714 };
73  static float p3[5] = { (float).1788544503,(float).09359161662,(float)
74  .006325387654,(float)6.611667319e-5,(float)-2.031049101e-6 };
75  static float q3[5] = { (float)1.,(float).6097809921,(float).2560616665,(
76  float).04746722384,(float).006957301675 };
77  static float p4[5] = { (float).9874054407,(float)118.6723273,(float)
78  849.279436,(float)-743.7792444,(float)427.0262186 };
79  static float q4[5] = { (float)1.,(float)106.8615961,(float)337.6496214,(
80  float)2016.712389,(float)1597.063511 };
81  static float p5[5] = { (float)1.003675074,(float)167.5702434,(float)
82  4789.711289,(float)21217.86767,(float)-22324.9491 };
83 
84  /* System generated locals */
85  float ret_val, r__1;
86 
87  /* Local variables */
88  static float u, v;
89  v = x;
90  if (v < (float)-5.5) {
91  u = exp(v + (float)1.);
92  ret_val = exp(-1 / u) / sqrt(u) * (float).3989422803 * ((a1[0] + (a1[
93  1] + a1[2] * u) * u) * u + 1);
94  } else if (v < (float)-1.) {
95  u = exp(-v - 1);
96  ret_val = exp(-u) * sqrt(u) * (p1[0] + (p1[1] + (p1[2] + (p1[3] + p1[
97  4] * v) * v) * v) * v) / (q1[0] + (q1[1] + (q1[2] + (q1[3] +
98  q1[4] * v) * v) * v) * v);
99  } else if (v < (float)1.) {
100  ret_val = (p2[0] + (p2[1] + (p2[2] + (p2[3] + p2[4] * v) * v) * v) *
101  v) / (q2[0] + (q2[1] + (q2[2] + (q2[3] + q2[4] * v) * v) * v)
102  * v);
103  } else if (v < (float)5.) {
104  ret_val = (p3[0] + (p3[1] + (p3[2] + (p3[3] + p3[4] * v) * v) * v) *
105  v) / (q3[0] + (q3[1] + (q3[2] + (q3[3] + q3[4] * v) * v) * v)
106  * v);
107  } else if (v < (float)12.) {
108  u = 1 / v;
109 /* Computing 2nd power */
110  r__1 = u;
111  ret_val = r__1 * r__1 * (p4[0] + (p4[1] + (p4[2] + (p4[3] + p4[4] * u)
112  * u) * u) * u) / (q4[0] + (q4[1] + (q4[2] + (q4[3] + q4[4] *
113  u) * u) * u) * u);
114  } else if (v < (float)50.) {
115  u = 1 / v;
116 /* Computing 2nd power */
117  r__1 = u;
118  ret_val = r__1 * r__1 * (p5[0] + (p5[1] + (p5[2] + (p5[3] + p5[4] * u)
119  * u) * u) * u) / (q5[0] + (q5[1] + (q5[2] + (q5[3] + q5[4] *
120  u) * u) * u) * u);
121  } else if (v < (float)300.) {
122  u = 1 / v;
123 /* Computing 2nd power */
124  r__1 = u;
125  ret_val = r__1 * r__1 * (p6[0] + (p6[1] + (p6[2] + (p6[3] + p6[4] * u)
126  * u) * u) * u) / (q6[0] + (q6[1] + (q6[2] + (q6[3] + q6[4] *
127  u) * u) * u) * u);
128  } else {
129  u = 1 / (v - v * log(v) / (v + 1));
130 /* Computing 2nd power */
131  r__1 = u;
132  ret_val = r__1 * r__1 * ((a2[0] + a2[1] * u) * u + 1);
133  }
134  return ret_val;
135 
136 }
137 
138 } // namespace Genfun
Genfun::Parameter::getValue
virtual double getValue() const
Definition: Parameter.cc:27
Genfun::Landau
Definition: CLHEP/GenericFunctions/Landau.hh:20
Genfun::AbsFunction
Definition: CLHEP/GenericFunctions/AbsFunction.hh:48
Landau.hh
Variable.hh
v
they are gone ZOOM Features Discontinued The following features of the ZOOM package were felt to be extreme overkill These have been after checking that no existing user code was utilizing as in SpaceVector v
Definition: keyMergeIssues.doc:324
Genfun::Landau::width
Parameter & width()
Definition: Landau.cc:41
s
Methods applicble to containers of as in std::list< LorentzVector > s
Definition: keyMergeIssues.doc:328
Genfun::Landau::operator()
virtual double operator()(double argument) const
Definition: Landau.cc:30
Genfun::Parameter
Definition: CLHEP/GenericFunctions/Parameter.hh:35
x
any side effects of that construction would occur twice The semantics of throw x
Definition: whyZMthrowRethrows.txt:37
Genfun::Landau::peak
Parameter & peak()
Definition: Landau.cc:37
FUNCTION_OBJECT_IMP
#define FUNCTION_OBJECT_IMP(classname)
Definition: CLHEP/GenericFunctions/AbsFunction.hh:156
Genfun
Definition: CLHEP/GenericFunctions/Abs.hh:14