]> git.uio.no Git - u/mrichter/AliRoot.git/blob - spdf/sasbeh.F
Update master to aliroot
[u/mrichter/AliRoot.git] / spdf / sasbeh.F
1 #include "pdf/pilot.h"
2       SUBROUTINE SASBEH(KF,X,Q2,P2,PM2,XPBH)
3 C...Purpose: to evaluate the Bethe-Heitler cross section for
4 C...heavy flavour production.
5       DATA AEM2PI/0.0011614/
6  
7 C...Reset output.
8       XPBH=0.
9       SIGBH=0.
10  
11 C...Check kinematics limits.
12       IF(X.GE.Q2/(4.*PM2+Q2+P2)) RETURN
13       W2=Q2*(1.-X)/X-P2
14       BETA2=1.-4.*PM2/W2
15       IF(BETA2.LT.1E-10) RETURN
16       RMQ=4.*PM2/Q2
17  
18 C...Simple case: P2 = 0.
19       IF(P2.LT.1E-4) THEN
20         BETA=SQRT(BETA2)
21         IF(BETA.LT.0.99) THEN
22           XBL=LOG((1.+BETA)/(1.-BETA))
23         ELSE
24           XBL=LOG((1.+BETA)**2*W2/(4.*PM2))
25         ENDIF
26         SIGBH=BETA*(8.*X*(1.-X)-1.-RMQ*X*(1.-X))+
27      &  XBL*(X**2+(1.-X)**2+RMQ*X*(1.-3.*X)-0.5*RMQ**2*X**2)
28  
29 C...Complicated case: P2 > 0, based on approximation of
30 C...C.T. Hill and G.G. Ross, Nucl. Phys. B148 (1979) 373
31       ELSE
32         RPQ=1.-4.*X**2*P2/Q2
33         IF(RPQ.GT.1E-10) THEN
34           RPBE=SQRT(RPQ*BETA2)
35           IF(RPBE.LT.0.99) THEN
36             XBL=LOG((1.+RPBE)/(1.-RPBE))
37             XBI=2.*RPBE/(1.-RPBE**2)
38           ELSE
39             RPBESN=4.*PM2/W2+(4.*X**2*P2/Q2)*BETA2
40             XBL=LOG((1.+RPBE)**2/RPBESN)
41             XBI=2.*RPBE/RPBESN
42           ENDIF
43           SIGBH=BETA*(6.*X*(1.-X)-1.)+
44      &    XBL*(X**2+(1.-X)**2+RMQ*X*(1.-3.*X)-0.5*RMQ**2*X**2)+
45      &    XBI*(2.*X/Q2)*(PM2*X*(2.-RMQ)-P2*X)
46         ENDIF
47       ENDIF
48  
49 C...Multiply by charge-squared etc. to get parton distribution.
50       CHSQ=1./9.
51       IF(IABS(KF).EQ.2.OR.IABS(KF).EQ.4) CHSQ=4./9.
52       XPBH=3.*CHSQ*AEM2PI*X*SIGBH
53  
54       RETURN
55       END