]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PDF/spdf/sasgam2.F
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
Do not unload gAlice, it is needed until the end of the simulation run
[u/mrichter/AliRoot.git] / PDF / spdf / sasgam2.F
CommitLineData
21886bb6 1#include "pdf/pilot.h"
2C...SaSgam version 2 - parton distributions of the photon
3C...by Gerhard A. Schuler and Torbjorn Sjostrand
4C...For further information see Z. Phys. C68 (1995) 607
5C...and CERN-TH/96-04 and LU TP 96-2.
6C...Program last changed on 18 January 1996.
7
8C!!!Note that one further call parameter - IP2 - has been added
9C!!!to the SASGAM argument list compared with version 1.
10
11C...The user should only need to call the SASGAM routine,
12C...which in turn calls the auxiliary routines SASVMD, SASANO,
13C...SASBEH and SASDIR. The package is self-contained.
14
15C...One particular aspect of these parametrizations is that F2 for
16C...the photon is not obtained just as the charge-squared-weighted
17C...sum of quark distributions, but differ in the treatment of
18C...heavy flavours (in F2 the DIS relation W2 = Q2*(1-x)/x restricts
19C...the kinematics range of heavy-flavour production, but the same
20C...kinematics is not relevant e.g. for jet production) and, for the
21C...'MSbar' fits, in the addition of a Cgamma term related to the
22C...separation of direct processes. Schematically:
23C...PDF = VMD (rho, omega, phi) + anomalous (d, u, s, c, b).
24C...F2 = VMD (rho, omega, phi) + anomalous (d, u, s) +
25C... Bethe-Heitler (c, b) (+ Cgamma (d, u, s)).
26C...The J/psi and Upsilon states have not been included in the VMD sum,
27C...but low c and b masses in the other components should compensate
28C...for this in a duality sense.
29
30C...The calling sequence is the following:
31C CALL SASGAM2(ISET,X,Q2,P2,IP2,F2GM,XPDFGM)
32C...with the following declaration statement:
33C DIMENSION XPDFGM(-6:6)
34C...and, optionally, further information in:
35C COMMON/SASCOM/XPVMD(-6:6),XPANL(-6:6),XPANH(-6:6),XPBEH(-6:6),
36C &XPDIR(-6:6)
37C COMMON/SASVAL/VXPVMD(-6:6),VXPANL(-6:6),VXPANH(-6:6),VXPDGM(-6:6)
38C...Input: ISET = 1 : SaS set 1D ('DIS', Q0 = 0.6 GeV)
39C = 2 : SaS set 1M ('MSbar', Q0 = 0.6 GeV)
40C = 3 : SaS set 2D ('DIS', Q0 = 2 GeV)
41C = 4 : SaS set 2M ('MSbar', Q0 = 2 GeV)
42C X : x value.
43C Q2 : Q2 value.
44C P2 : P2 value; should be = 0. for an on-shell photon.
45C IP2 : scheme used to evaluate off-shell anomalous component.
46C = 0 : recommended default, see = 7.
47C = 1 : dipole dampening by integration; very time-consuming.
48C = 2 : P_0^2 = max( Q_0^2, P^2 )
49C = 3 : P'_0^2 = Q_0^2 + P^2.
50C = 4 : P_{eff} that preserves momentum sum.
51C = 5 : P_{int} that preserves momentum and average
52C evolution range.
53C = 6 : P_{eff}, matched to P_0 in P2 -> Q2 limit.
54C = 7 : P_{eff}, matched to P_0 in P2 -> Q2 limit.
55C...Output: F2GM : F2 value of the photon (including factors of alpha_em).
56C XPFDGM : x times parton distribution functions of the photon,
57C with elements 0 = g, 1 = d, 2 = u, 3 = s, 4 = c, 5 = b,
58C 6 = t (always empty!), - for antiquarks (result is same).
59C...The breakdown by component is stored in the commonblock SASCOM,
60C with elements as above.
61C XPVMD : rho, omega, phi VMD part only of output.
62C XPANL : d, u, s anomalous part only of output.
63C XPANH : c, b anomalous part only of output.
64C XPBEH : c, b Bethe-Heitler part only of output.
65C XPDIR : Cgamma (direct contribution) part only of output.
66C...The above arrays do not distinguish valence and sea contributions,
67C...although this information is available internally. The additional
68C...commonblock SASVAL provides the valence part only of the above
69C...distributions. Array names VXPVMD, VXPANL and VXPANH correspond
70C...to XPVMD, XPANL and XPANH, while XPBEH and XPDIR are valence only
71C...and therefore not given doubly. VXPDGM gives the sum of valence
72C...parts, and so matches XPDFGM. The difference, i.e. XPVMD-VXPVMD
73C...and so on, gives the sea part only.
74
75 SUBROUTINE SASGAM2(ISET,X,Q2,P2,IP2,F2GM,XPDFGM)
76C...Purpose: to construct the F2 and parton distributions of the photon
77C...by summing homogeneous (VMD) and inhomogeneous (anomalous) terms.
78C...For F2, c and b are included by the Bethe-Heitler formula;
79C...in the 'MSbar' scheme additionally a Cgamma term is added.
80 DIMENSION XPDFGM(-6:6)
81 COMMON/SASCOM/XPVMD(-6:6),XPANL(-6:6),XPANH(-6:6),XPBEH(-6:6),
82 &XPDIR(-6:6)
83 COMMON/SASVAL/VXPVMD(-6:6),VXPANL(-6:6),VXPANH(-6:6),VXPDGM(-6:6)
84 SAVE /SASCOM/,/SASVAL/
85
86C...Temporary array.
87 DIMENSION XPGA(-6:6), VXPGA(-6:6)
88C...Charm and bottom masses (low to compensate for J/psi etc.).
89 DATA PMC/1.3/, PMB/4.6/
90C...alpha_em and alpha_em/(2*pi).
91 DATA AEM/0.007297/, AEM2PI/0.0011614/
92C...Lambda value for 4 flavours.
93 DATA ALAM/0.20/
94C...Mixture u/(u+d), = 0.5 for incoherent and = 0.8 for coherent sum.
95 DATA FRACU/0.8/
96C...VMD couplings f_V**2/(4*pi).
97 DATA FRHO/2.20/, FOMEGA/23.6/, FPHI/18.4/
98C...Masses for rho (=omega) and phi.
99 DATA PMRHO/0.770/, PMPHI/1.020/
100C...Number of points in integration for IP2=1.
101 DATA NSTEP/100/
102
103C...Reset output.
104 F2GM=0.
105 DO 100 KFL=-6,6
106 XPDFGM(KFL)=0.
107 XPVMD(KFL)=0.
108 XPANL(KFL)=0.
109 XPANH(KFL)=0.
110 XPBEH(KFL)=0.
111 XPDIR(KFL)=0.
112 VXPVMD(KFL)=0.
113 VXPANL(KFL)=0.
114 VXPANH(KFL)=0.
115 VXPDGM(KFL)=0.
116 100 CONTINUE
117
118C...Check that input sensible.
119 IF(ISET.LE.0.OR.ISET.GE.5) THEN
120 WRITE(*,*) ' FATAL ERROR: SaSgam called for unknown set'
121 WRITE(*,*) ' ISET = ',ISET
122 STOP
123 ENDIF
124 IF(X.LE.0..OR.X.GT.1.) THEN
125 WRITE(*,*) ' FATAL ERROR: SaSgam called for unphysical x'
126 WRITE(*,*) ' X = ',X
127 STOP
128 ENDIF
129
130C...Set Q0 cut-off parameter as function of set used.
131 IF(ISET.LE.2) THEN
132 Q0=0.6
133 ELSE
134 Q0=2.
135 ENDIF
136 Q02=Q0**2
137
138C...Scale choice for off-shell photon; common factors.
139 Q2A=Q2
140 FACNOR=1.
141 IF(IP2.EQ.1) THEN
142 P2MX=P2+Q02
143 Q2A=Q2+P2*Q02/MAX(Q02,Q2)
144 FACNOR=LOG(Q2/Q02)/NSTEP
145 ELSEIF(IP2.EQ.2) THEN
146 P2MX=MAX(P2,Q02)
147 ELSEIF(IP2.EQ.3) THEN
148 P2MX=P2+Q02
149 Q2A=Q2+P2*Q02/MAX(Q02,Q2)
150 ELSEIF(IP2.EQ.4) THEN
151 P2MX=Q2*(Q02+P2)/(Q2+P2)*EXP(P2*(Q2-Q02)/
152 & ((Q2+P2)*(Q02+P2)))
153 ELSEIF(IP2.EQ.5) THEN
154 P2MXA=Q2*(Q02+P2)/(Q2+P2)*EXP(P2*(Q2-Q02)/
155 & ((Q2+P2)*(Q02+P2)))
156 P2MX=Q0*SQRT(P2MXA)
157 FACNOR=LOG(Q2/P2MXA)/LOG(Q2/P2MX)
158 ELSEIF(IP2.EQ.6) THEN
159 P2MX=Q2*(Q02+P2)/(Q2+P2)*EXP(P2*(Q2-Q02)/
160 & ((Q2+P2)*(Q02+P2)))
161 P2MX=MAX(0.,1.-P2/Q2)*P2MX+MIN(1.,P2/Q2)*MAX(P2,Q02)
162 ELSE
163 P2MXA=Q2*(Q02+P2)/(Q2+P2)*EXP(P2*(Q2-Q02)/
164 & ((Q2+P2)*(Q02+P2)))
165 P2MX=Q0*SQRT(P2MXA)
166 P2MXB=P2MX
167 P2MX=MAX(0.,1.-P2/Q2)*P2MX+MIN(1.,P2/Q2)*MAX(P2,Q02)
168 P2MXB=MAX(0.,1.-P2/Q2)*P2MXB+MIN(1.,P2/Q2)*P2MXA
169 FACNOR=LOG(Q2/P2MXA)/LOG(Q2/P2MXB)
170 ENDIF
171
172C...Call VMD parametrization for d quark and use to give rho, omega,
173C...phi. Note dipole dampening for off-shell photon.
174 CALL SASVMD(ISET,1,X,Q2A,P2MX,ALAM,XPGA,VXPGA)
175 XFVAL=VXPGA(1)
176 XPGA(1)=XPGA(2)
177 XPGA(-1)=XPGA(-2)
178 FACUD=AEM*(1./FRHO+1./FOMEGA)*(PMRHO**2/(PMRHO**2+P2))**2
179 FACS=AEM*(1./FPHI)*(PMPHI**2/(PMPHI**2+P2))**2
180 DO 110 KFL=-5,5
181 XPVMD(KFL)=(FACUD+FACS)*XPGA(KFL)
182 110 CONTINUE
183 XPVMD(1)=XPVMD(1)+(1.-FRACU)*FACUD*XFVAL
184 XPVMD(2)=XPVMD(2)+FRACU*FACUD*XFVAL
185 XPVMD(3)=XPVMD(3)+FACS*XFVAL
186 XPVMD(-1)=XPVMD(-1)+(1.-FRACU)*FACUD*XFVAL
187 XPVMD(-2)=XPVMD(-2)+FRACU*FACUD*XFVAL
188 XPVMD(-3)=XPVMD(-3)+FACS*XFVAL
189 VXPVMD(1)=(1.-FRACU)*FACUD*XFVAL
190 VXPVMD(2)=FRACU*FACUD*XFVAL
191 VXPVMD(3)=FACS*XFVAL
192 VXPVMD(-1)=(1.-FRACU)*FACUD*XFVAL
193 VXPVMD(-2)=FRACU*FACUD*XFVAL
194 VXPVMD(-3)=FACS*XFVAL
195
196 IF(IP2.NE.1) THEN
197C...Anomalous parametrizations for different strategies
198C...for off-shell photons; except full integration.
199
200C...Call anomalous parametrization for d + u + s.
201 CALL SASANO(-3,X,Q2A,P2MX,ALAM,XPGA,VXPGA)
202 DO 120 KFL=-5,5
203 XPANL(KFL)=FACNOR*XPGA(KFL)
204 VXPANL(KFL)=FACNOR*VXPGA(KFL)
205 120 CONTINUE
206
207C...Call anomalous parametrization for c and b.
208 CALL SASANO(4,X,Q2A,P2MX,ALAM,XPGA,VXPGA)
209 DO 130 KFL=-5,5
210 XPANH(KFL)=FACNOR*XPGA(KFL)
211 VXPANH(KFL)=FACNOR*VXPGA(KFL)
212 130 CONTINUE
213 CALL SASANO(5,X,Q2A,P2MX,ALAM,XPGA,VXPGA)
214 DO 140 KFL=-5,5
215 XPANH(KFL)=XPANH(KFL)+FACNOR*XPGA(KFL)
216 VXPANH(KFL)=VXPANH(KFL)+FACNOR*VXPGA(KFL)
217 140 CONTINUE
218
219 ELSE
220C...Special option: loop over flavours and integrate over k2.
221 DO 170 KF=1,5
222 DO 160 ISTEP=1,NSTEP
223 Q2STEP=Q02*(Q2/Q02)**((ISTEP-0.5)/NSTEP)
224 IF((KF.EQ.4.AND.Q2STEP.LT.PMC**2).OR.
225 & (KF.EQ.5.AND.Q2STEP.LT.PMB**2)) GOTO 160
226 CALL SASVMD(0,KF,X,Q2,Q2STEP,ALAM,XPGA,VXPGA)
227 FACQ=AEM2PI*(Q2STEP/(Q2STEP+P2))**2*FACNOR
228 IF(MOD(KF,2).EQ.0) FACQ=FACQ*(8./9.)
229 IF(MOD(KF,2).EQ.1) FACQ=FACQ*(2./9.)
230 DO 150 KFL=-5,5
231 IF(KF.LE.3) XPANL(KFL)=XPANL(KFL)+FACQ*XPGA(KFL)
232 IF(KF.GE.4) XPANH(KFL)=XPANH(KFL)+FACQ*XPGA(KFL)
233 IF(KF.LE.3) VXPANL(KFL)=VXPANL(KFL)+FACQ*VXPGA(KFL)
234 IF(KF.GE.4) VXPANH(KFL)=VXPANH(KFL)+FACQ*VXPGA(KFL)
235 150 CONTINUE
236 160 CONTINUE
237 170 CONTINUE
238 ENDIF
239
240C...Call Bethe-Heitler term expression for charm and bottom.
241 CALL SASBEH(4,X,Q2,P2,PMC**2,XPBH)
242 XPBEH(4)=XPBH
243 XPBEH(-4)=XPBH
244 CALL SASBEH(5,X,Q2,P2,PMB**2,XPBH)
245 XPBEH(5)=XPBH
246 XPBEH(-5)=XPBH
247
248C...For MSbar subtraction call C^gamma term expression for d, u, s.
249 IF(ISET.EQ.2.OR.ISET.EQ.4) THEN
250 CALL SASDIR(X,Q2,P2,Q02,XPGA)
251 DO 180 KFL=-5,5
252 XPDIR(KFL)=XPGA(KFL)
253 180 CONTINUE
254 ENDIF
255
256C...Store result in output array.
257 DO 190 KFL=-5,5
258 CHSQ=1./9.
259 IF(IABS(KFL).EQ.2.OR.IABS(KFL).EQ.4) CHSQ=4./9.
260 XPF2=XPVMD(KFL)+XPANL(KFL)+XPBEH(KFL)+XPDIR(KFL)
261 IF(KFL.NE.0) F2GM=F2GM+CHSQ*XPF2
262 XPDFGM(KFL)=XPVMD(KFL)+XPANL(KFL)+XPANH(KFL)
263 VXPDGM(KFL)=VXPVMD(KFL)+VXPANL(KFL)+VXPANH(KFL)
264 190 CONTINUE
265
266 RETURN
267 END