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
10using namespace std;
11
12namespace Genfun {
14
16 _peak("Peak", 5.0 ,0, 10),
17 _width("Width",1.0,0, 10)
18{}
19
21}
22
23Landau::Landau(const Landau & right):
24AbsFunction(right),
25_peak(right._peak),
26_width(right._width)
27{
28}
29
30double 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
45const Parameter & Landau::peak() const {
46 return _peak;
47}
48
49const Parameter & Landau::width() const {
50 return _width;
51}
52
53double 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
#define FUNCTION_OBJECT_IMP(classname)
virtual ~Landau()
Definition: Landau.cc:20
Parameter & peak()
Definition: Landau.cc:37
virtual double operator()(double argument) const
Definition: Landau.cc:30
Parameter & width()
Definition: Landau.cc:41
virtual double getValue() const
Definition: Parameter.cc:27