]>
Commit | Line | Data |
---|---|---|
0ca57c2f | 1 | // $Id: EvtDiLog.cpp,v 1.3 2009-03-16 15:52:42 robbep Exp $ |
da0e9ce3 | 2 | // Include files |
3 | ||
4 | // local | |
5 | #include "EvtGenBase/EvtDiLog.hh" | |
6 | ||
7 | #include <cmath> | |
8 | ||
9 | //----------------------------------------------------------------------------- | |
10 | // Implementation file for class : EvtDiLog | |
11 | // | |
12 | // 2007-01-23 : Patrick Robbe | |
13 | //----------------------------------------------------------------------------- | |
14 | ||
15 | double EvtDiLog::DiLog( double x ) { | |
16 | ||
17 | double h , t , y , s , a , alfa , b0, b1, b2 ; | |
18 | if ( x == 1. ) h = PI6 ; | |
19 | else if ( x == -1. ) h = -PI12 ; | |
20 | else { | |
21 | t = -x ; | |
22 | if ( t <= -2. ) { | |
23 | y = -1./(1.+t) ; | |
24 | s = 1. ; | |
25 | a = -PI3 + HF * ( std::pow( log(-t) , 2 ) - | |
26 | std::pow( log( 1. + 1./t ) , 2 ) ) ; | |
27 | } else if ( t < -1. ) { | |
28 | y = -1. - t ; | |
29 | s = -1. ; | |
30 | a = log( -t ) ; | |
31 | a = -PI6 + a * ( a + log( 1. + 1./t ) ) ; | |
32 | } else if ( t <= -HF ) { | |
33 | y = - (1. + t ) / t ; | |
34 | s = 1. ; | |
35 | a = log( -t ) ; | |
36 | a = -PI6 + a * ( -HF * a + log( 1. + t ) ) ; | |
37 | } else if ( t < 0 ) { | |
38 | y = -t / ( 1. + t ) ; | |
39 | s = -1. ; | |
40 | a = HF * std::pow( log( 1. + t ) , 2 ) ; | |
41 | } else if ( t <= 1. ) { | |
42 | y = t ; | |
43 | s = 1. ; | |
44 | a = 0. ; | |
45 | } else { | |
46 | y = 1. / t ; | |
47 | s = -1. ; | |
48 | a = PI6 + HF * std::pow( log( t ) , 2 ) ; | |
49 | } | |
50 | ||
51 | h = y + y - 1. ; | |
52 | alfa = h + h ; | |
53 | b1 = 0. ; | |
54 | b2 = 0. ; | |
55 | for ( int i = 19 ; i >= 0 ; --i ) { | |
56 | b0 = C[ i ] + alfa * b1 - b2 ; | |
57 | b2 = b1 ; | |
58 | b1 = b0 ; | |
59 | } | |
60 | ||
61 | h = -(s * ( b0 -h * b2 ) + a ) ; | |
62 | } | |
63 | ||
64 | return h ; | |
65 | } |