forked from sPHENIX-Collaboration/coresoftware
/
pyhispen.f
executable file
·86 lines (75 loc) · 2.49 KB
/
pyhispen.f
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
C***********************************************************************
FUNCTION PYHISPEN(XREIN,XIMIN,IREIM)
C...Calculates real and imaginary part of Spence function; see
C...G. 't Hooft and M. Veltman, Nucl. Phys. B153 (1979) 365.
COMMON/LUDAT1/MSTU(200),PARU(200),MSTJ(200),PARJ(200)
SAVE /LUDAT1/
DIMENSION B(0:14)
DATA B/
& 1.000000E+00, -5.000000E-01, 1.666667E-01,
& 0.000000E+00, -3.333333E-02, 0.000000E+00,
& 2.380952E-02, 0.000000E+00, -3.333333E-02,
& 0.000000E+00, 7.575757E-02, 0.000000E+00,
&-2.531135E-01, 0.000000E+00, 1.166667E+00/
XRE=XREIN
XIM=XIMIN
IF(ABS(1.-XRE).LT.1.E-6.AND.ABS(XIM).LT.1.E-6) THEN
IF(IREIM.EQ.1) PYHISPEN=PARU(1)**2/6.
IF(IREIM.EQ.2) PYHISPEN=0.
RETURN
ENDIF
XMOD=SQRT(XRE**2+XIM**2)
IF(XMOD.LT.1.E-6) THEN
IF(IREIM.EQ.1) PYHISPEN=0.
IF(IREIM.EQ.2) PYHISPEN=0.
RETURN
ENDIF
XARG=SIGN(ACOS(XRE/XMOD),XIM)
SP0RE=0.
SP0IM=0.
SGN=1.
IF(XMOD.GT.1.) THEN
ALGXRE=LOG(XMOD)
ALGXIM=XARG-SIGN(PARU(1),XARG)
SP0RE=-PARU(1)**2/6.-(ALGXRE**2-ALGXIM**2)/2.
SP0IM=-ALGXRE*ALGXIM
SGN=-1.
XMOD=1./XMOD
XARG=-XARG
XRE=XMOD*COS(XARG)
XIM=XMOD*SIN(XARG)
ENDIF
IF(XRE.GT.0.5) THEN
ALGXRE=LOG(XMOD)
ALGXIM=XARG
XRE=1.-XRE
XIM=-XIM
XMOD=SQRT(XRE**2+XIM**2)
XARG=SIGN(ACOS(XRE/XMOD),XIM)
ALGYRE=LOG(XMOD)
ALGYIM=XARG
SP0RE=SP0RE+SGN*(PARU(1)**2/6.-(ALGXRE*ALGYRE-ALGXIM*ALGYIM))
SP0IM=SP0IM-SGN*(ALGXRE*ALGYIM+ALGXIM*ALGYRE)
SGN=-SGN
ENDIF
XRE=1.-XRE
XIM=-XIM
XMOD=SQRT(XRE**2+XIM**2)
XARG=SIGN(ACOS(XRE/XMOD),XIM)
ZRE=-LOG(XMOD)
ZIM=-XARG
SPRE=0.
SPIM=0.
SAVERE=1.
SAVEIM=0.
DO 100 I=0,14
TERMRE=(SAVERE*ZRE-SAVEIM*ZIM)/FLOAT(I+1)
TERMIM=(SAVERE*ZIM+SAVEIM*ZRE)/FLOAT(I+1)
SAVERE=TERMRE
SAVEIM=TERMIM
SPRE=SPRE+B(I)*TERMRE
100 SPIM=SPIM+B(I)*TERMIM
IF(IREIM.EQ.1) PYHISPEN=SP0RE+SGN*SPRE
IF(IREIM.EQ.2) PYHISPEN=SP0IM+SGN*SPIM
RETURN
END