[u/mrichter/AliRoot.git] / HERWIG / src / sasgam.f


CDECK  ID>, STRUCTM.

*CMZ :S        E26/04/91  11.11.54  by  Bryan Webber

*-- Author :    Bryan Webber

C-----------------------------------------------------------------------

C      SUBROUTINE STRUCTM(X,QSCA,UPV,DNV,USEA,DSEA,STR,CHM,BOT,TOP,GLU)

C-----------------------------------------------------------------------

C     DUMMY SUBROUTINE: DELETE IF YOU USE PDFLIB CERN-LIBRARY

C     PACKAGE FOR NUCLEON STRUCTURE FUNCTIONS

C-----------------------------------------------------------------------

C      DOUBLE PRECISION X,QSCA,UPV,DNV,USEA,DSEA,STR,CHM,BOT,TOP,GLU

C      WRITE (6,10)

C  10  FORMAT(/10X,'STRUCTM CALLED BUT NOT LINKED')

C      STOP

C      END

C-----------------------------------------------------------------------

C...SaSgam version 2 - parton distributions of the photon

C...by Gerhard A. Schuler and Torbjorn Sjostrand

C...For further information see Z. Phys. C68 (1995) 607

C...and CERN-TH/96-04 and LU TP 96-2.

C...Program last changed on 18 January 1996.

C

C!!!Note that one further call parameter - IP2 - has been added

C!!!to the SASGAM argument list compared with version 1.

C

C...The user should only need to call the SASGAM routine,

C...which in turn calls the auxiliary routines SASVMD, SASANO,

C...SASBEH and SASDIR. The package is self-contained.

C

C...One particular aspect of these parametrizations is that F2 for

C...the photon is not obtained just as the charge-squared-weighted

C...sum of quark distributions, but differ in the treatment of

C...heavy flavours (in F2 the DIS relation W2 = Q2*(1-x)/x restricts

C...the kinematics range of heavy-flavour production, but the same

C...kinematics is not relevant e.g. for jet production) and, for the

C...'MSbar' fits, in the addition of a Cgamma term related to the

C...separation of direct processes. Schematically:

C...PDF = VMD (rho, omega, phi) + anomalous (d, u, s, c, b).

C...F2  = VMD (rho, omega, phi) + anomalous (d, u, s) +

C...      Bethe-Heitler (c, b) (+ Cgamma (d, u, s)).

C...The J/psi and Upsilon states have not been included in the VMD sum,

C...but low c and b masses in the other components should compensate

C...for this in a duality sense.

C

C...The calling sequence is the following:

C     CALL SASGAM(ISET,X,Q2,P2,IP2,F2GM,XPDFGM)

C...with the following declaration statement:

C     DIMENSION XPDFGM(-6:6)

C...and, optionally, further information in:

C     COMMON/SASCOM/XPVMD(-6:6),XPANL(-6:6),XPANH(-6:6),XPBEH(-6:6),

C    &XPDIR(-6:6)

C     COMMON/SASVAL/VXPVMD(-6:6),VXPANL(-6:6),VXPANH(-6:6),VXPDGM(-6:6)

C...Input:  ISET = 1 : SaS set 1D ('DIS',   Q0 = 0.6 GeV)

C                = 2 : SaS set 1M ('MSbar', Q0 = 0.6 GeV)

C                = 3 : SaS set 2D ('DIS',   Q0 =  2  GeV)

C                = 4 : SaS set 2M ('MSbar', Q0 =  2  GeV)

C           X : x value.

C           Q2 : Q2 value.

C           P2 : P2 value; should be = 0. for an on-shell photon.

C           IP2 : scheme used to evaluate off-shell anomalous component.

C               = 0 : recommended default, see = 7.

C               = 1 : dipole dampening by integration; very time-consuming.

C               = 2 : P_0^2 = max( Q_0^2, P^2 )

C               = 3 : P'_0^2 = Q_0^2 + P^2.

C               = 4 : P_{eff} that preserves momentum sum.

C               = 5 : P_{int} that preserves momentum and average

C                     evolution range.

C               = 6 : P_{eff}, matched to P_0 in P2 -> Q2 limit.

C               = 7 : P_{eff}, matched to P_0 in P2 -> Q2 limit.

C...Output: F2GM : F2 value of the photon (including factors of alpha_em).

C           XPFDGM :  x times parton distribution functions of the photon,

C               with elements 0 = g, 1 = d, 2 = u, 3 = s, 4 = c, 5 = b,

C               6 = t (always empty!), - for antiquarks (result is same).

C...The breakdown by component is stored in the commonblock SASCOM,

C               with elements as above.

C           XPVMD : rho, omega, phi VMD part only of output.

C           XPANL : d, u, s anomalous part only of output.

C           XPANH : c, b anomalous part only of output.

C           XPBEH : c, b Bethe-Heitler part only of output.

C           XPDIR : Cgamma (direct contribution) part only of output.

C...The above arrays do not distinguish valence and sea contributions,

C...although this information is available internally. The additional

C...commonblock SASVAL provides the valence part only of the above

C...distributions. Array names VXPVMD, VXPANL and VXPANH correspond

C...to XPVMD, XPANL and XPANH, while XPBEH and XPDIR are valence only

C...and therefore not given doubly. VXPDGM gives the sum of valence

C...parts, and so matches XPDFGM. The difference, i.e. XPVMD-VXPVMD

C...and so on, gives the sea part only.

C

      SUBROUTINE SASGAM(ISET,X,Q2,P2,IP2,F2GM,XPDFGM)

C...Purpose: to construct the F2 and parton distributions of the photon

C...by summing homogeneous (VMD) and inhomogeneous (anomalous) terms.

C...For F2, c and b are included by the Bethe-Heitler formula;

C...in the 'MSbar' scheme additionally a Cgamma term is added.

      DIMENSION XPDFGM(-6:6)

      COMMON/SASCOM/XPVMD(-6:6),XPANL(-6:6),XPANH(-6:6),XPBEH(-6:6),

     &XPDIR(-6:6)

      COMMON/SASVAL/VXPVMD(-6:6),VXPANL(-6:6),VXPANH(-6:6),VXPDGM(-6:6)

      SAVE /SASCOM/,/SASVAL/

C

C...Temporary array.

      DIMENSION XPGA(-6:6), VXPGA(-6:6)

C...Charm and bottom masses (low to compensate for J/psi etc.).

      DATA PMC/1.3/, PMB/4.6/

C...alpha_em and alpha_em/(2*pi).

      DATA AEM/0.007297/, AEM2PI/0.0011614/

C...Lambda value for 4 flavours.

      DATA ALAM/0.20/

C...Mixture u/(u+d), = 0.5 for incoherent and = 0.8 for coherent sum.

      DATA FRACU/0.8/

C...VMD couplings f_V**2/(4*pi).

      DATA FRHO/2.20/, FOMEGA/23.6/, FPHI/18.4/

C...Masses for rho (=omega) and phi.

      DATA PMRHO/0.770/, PMPHI/1.020/

C...Number of points in integration for IP2=1.

      DATA NSTEP/100/

C

C...Reset output.

      F2GM=0.

      DO 100 KFL=-6,6

      XPDFGM(KFL)=0.

      XPVMD(KFL)=0.

      XPANL(KFL)=0.

      XPANH(KFL)=0.

      XPBEH(KFL)=0.

      XPDIR(KFL)=0.

      VXPVMD(KFL)=0.

      VXPANL(KFL)=0.

      VXPANH(KFL)=0.

      VXPDGM(KFL)=0.

  100 CONTINUE

C

C...Check that input sensible.

      IF(ISET.LE.0.OR.ISET.GE.5) THEN

        WRITE(*,*) ' FATAL ERROR: SaSgam called for unknown set'

        WRITE(*,*) ' ISET = ',ISET

        STOP

      ENDIF

      IF(X.LE.0..OR.X.GT.1.) THEN

        WRITE(*,*) ' FATAL ERROR: SaSgam called for unphysical x'

        WRITE(*,*) ' X = ',X

        STOP

      ENDIF

C

C...Set Q0 cut-off parameter as function of set used.

      IF(ISET.LE.2) THEN

        Q0=0.6

      ELSE

        Q0=2.

      ENDIF

      Q02=Q0**2

C

C...Scale choice for off-shell photon; common factors.

      Q2A=Q2

      FACNOR=1.

      IF(IP2.EQ.1) THEN

        P2MX=P2+Q02

        Q2A=Q2+P2*Q02/MAX(Q02,Q2)

        FACNOR=LOG(Q2/Q02)/NSTEP

      ELSEIF(IP2.EQ.2) THEN

        P2MX=MAX(P2,Q02)

      ELSEIF(IP2.EQ.3) THEN

        P2MX=P2+Q02

        Q2A=Q2+P2*Q02/MAX(Q02,Q2)

      ELSEIF(IP2.EQ.4) THEN

        P2MX=Q2*(Q02+P2)/(Q2+P2)*EXP(P2*(Q2-Q02)/

     &  ((Q2+P2)*(Q02+P2)))

      ELSEIF(IP2.EQ.5) THEN

        P2MXA=Q2*(Q02+P2)/(Q2+P2)*EXP(P2*(Q2-Q02)/

     &  ((Q2+P2)*(Q02+P2)))

        P2MX=Q0*SQRT(P2MXA)

        FACNOR=LOG(Q2/P2MXA)/LOG(Q2/P2MX)

      ELSEIF(IP2.EQ.6) THEN

        P2MX=Q2*(Q02+P2)/(Q2+P2)*EXP(P2*(Q2-Q02)/

     &  ((Q2+P2)*(Q02+P2)))

        P2MX=MAX(0.,1.-P2/Q2)*P2MX+MIN(1.,P2/Q2)*MAX(P2,Q02)

      ELSE

        P2MXA=Q2*(Q02+P2)/(Q2+P2)*EXP(P2*(Q2-Q02)/

     &  ((Q2+P2)*(Q02+P2)))

        P2MX=Q0*SQRT(P2MXA)

        P2MXB=P2MX

        P2MX=MAX(0.,1.-P2/Q2)*P2MX+MIN(1.,P2/Q2)*MAX(P2,Q02)

        P2MXB=MAX(0.,1.-P2/Q2)*P2MXB+MIN(1.,P2/Q2)*P2MXA

        FACNOR=LOG(Q2/P2MXA)/LOG(Q2/P2MXB)

      ENDIF

C

C...Call VMD parametrization for d quark and use to give rho, omega,

C...phi. Note dipole dampening for off-shell photon.

      CALL SASVMD(ISET,1,X,Q2A,P2MX,ALAM,XPGA,VXPGA)

      XFVAL=VXPGA(1)

      XPGA(1)=XPGA(2)

      XPGA(-1)=XPGA(-2)

      FACUD=AEM*(1./FRHO+1./FOMEGA)*(PMRHO**2/(PMRHO**2+P2))**2

      FACS=AEM*(1./FPHI)*(PMPHI**2/(PMPHI**2+P2))**2

      DO 110 KFL=-5,5

      XPVMD(KFL)=(FACUD+FACS)*XPGA(KFL)

  110 CONTINUE

      XPVMD(1)=XPVMD(1)+(1.-FRACU)*FACUD*XFVAL

      XPVMD(2)=XPVMD(2)+FRACU*FACUD*XFVAL

      XPVMD(3)=XPVMD(3)+FACS*XFVAL

      XPVMD(-1)=XPVMD(-1)+(1.-FRACU)*FACUD*XFVAL

      XPVMD(-2)=XPVMD(-2)+FRACU*FACUD*XFVAL

      XPVMD(-3)=XPVMD(-3)+FACS*XFVAL

      VXPVMD(1)=(1.-FRACU)*FACUD*XFVAL

      VXPVMD(2)=FRACU*FACUD*XFVAL

      VXPVMD(3)=FACS*XFVAL

      VXPVMD(-1)=(1.-FRACU)*FACUD*XFVAL

      VXPVMD(-2)=FRACU*FACUD*XFVAL

      VXPVMD(-3)=FACS*XFVAL

C

      IF(IP2.NE.1) THEN

C...Anomalous parametrizations for different strategies

C...for off-shell photons; except full integration.

C

C...Call anomalous parametrization for d + u + s.

        CALL SASANO(-3,X,Q2A,P2MX,ALAM,XPGA,VXPGA)

        DO 120 KFL=-5,5

        XPANL(KFL)=FACNOR*XPGA(KFL)

        VXPANL(KFL)=FACNOR*VXPGA(KFL)

  120   CONTINUE

C

C...Call anomalous parametrization for c and b.

        CALL SASANO(4,X,Q2A,P2MX,ALAM,XPGA,VXPGA)

        DO 130 KFL=-5,5

        XPANH(KFL)=FACNOR*XPGA(KFL)

        VXPANH(KFL)=FACNOR*VXPGA(KFL)

  130   CONTINUE

        CALL SASANO(5,X,Q2A,P2MX,ALAM,XPGA,VXPGA)

        DO 140 KFL=-5,5

        XPANH(KFL)=XPANH(KFL)+FACNOR*XPGA(KFL)

        VXPANH(KFL)=VXPANH(KFL)+FACNOR*VXPGA(KFL)

  140   CONTINUE

C

      ELSE

C...Special option: loop over flavours and integrate over k2.

        DO 170 KF=1,5

        DO 160 ISTEP=1,NSTEP

        Q2STEP=Q02*(Q2/Q02)**((ISTEP-0.5)/NSTEP)

        IF((KF.EQ.4.AND.Q2STEP.LT.PMC**2).OR.

     &  (KF.EQ.5.AND.Q2STEP.LT.PMB**2)) GOTO 160

        CALL SASVMD(0,KF,X,Q2,Q2STEP,ALAM,XPGA,VXPGA)

        FACQ=AEM2PI*(Q2STEP/(Q2STEP+P2))**2*FACNOR

        IF(MOD(KF,2).EQ.0) FACQ=FACQ*(8./9.)

        IF(MOD(KF,2).EQ.1) FACQ=FACQ*(2./9.)

        DO 150 KFL=-5,5

        IF(KF.LE.3) XPANL(KFL)=XPANL(KFL)+FACQ*XPGA(KFL)

        IF(KF.GE.4) XPANH(KFL)=XPANH(KFL)+FACQ*XPGA(KFL)

        IF(KF.LE.3) VXPANL(KFL)=VXPANL(KFL)+FACQ*VXPGA(KFL)

        IF(KF.GE.4) VXPANH(KFL)=VXPANH(KFL)+FACQ*VXPGA(KFL)

  150   CONTINUE

  160   CONTINUE

  170   CONTINUE

      ENDIF

C

C...Call Bethe-Heitler term expression for charm and bottom.

      CALL SASBEH(4,X,Q2,P2,PMC**2,XPBH)

      XPBEH(4)=XPBH

      XPBEH(-4)=XPBH

      CALL SASBEH(5,X,Q2,P2,PMB**2,XPBH)

      XPBEH(5)=XPBH

      XPBEH(-5)=XPBH

C

C...For MSbar subtraction call C^gamma term expression for d, u, s.

      IF(ISET.EQ.2.OR.ISET.EQ.4) THEN

        CALL SASDIR(X,Q2,P2,Q02,XPGA)

        DO 180 KFL=-5,5

        XPDIR(KFL)=XPGA(KFL)

  180   CONTINUE

      ENDIF

C

C...Store result in output array.

      DO 190 KFL=-5,5

      CHSQ=1./9.

      IF(IABS(KFL).EQ.2.OR.IABS(KFL).EQ.4) CHSQ=4./9.

      XPF2=XPVMD(KFL)+XPANL(KFL)+XPBEH(KFL)+XPDIR(KFL)

      IF(KFL.NE.0) F2GM=F2GM+CHSQ*XPF2

      XPDFGM(KFL)=XPVMD(KFL)+XPANL(KFL)+XPANH(KFL)

      VXPDGM(KFL)=VXPVMD(KFL)+VXPANL(KFL)+VXPANH(KFL)

  190 CONTINUE

C

      RETURN

      END
Commit	Line	Data
3820ca8e	1
	2	CDECK ID>, STRUCTM.
	3
	4	*CMZ :S E26/04/91 11.11.54 by Bryan Webber
	5
	6	*-- Author : Bryan Webber
	7
	8	C-----------------------------------------------------------------------
	9
	10	C SUBROUTINE STRUCTM(X,QSCA,UPV,DNV,USEA,DSEA,STR,CHM,BOT,TOP,GLU)
	11
	12	C-----------------------------------------------------------------------
	13
	14	C DUMMY SUBROUTINE: DELETE IF YOU USE PDFLIB CERN-LIBRARY
	15
	16	C PACKAGE FOR NUCLEON STRUCTURE FUNCTIONS
	17
	18	C-----------------------------------------------------------------------
	19
	20	C DOUBLE PRECISION X,QSCA,UPV,DNV,USEA,DSEA,STR,CHM,BOT,TOP,GLU
	21
	22	C WRITE (6,10)
	23
	24	C 10 FORMAT(/10X,'STRUCTM CALLED BUT NOT LINKED')
	25
	26	C STOP
	27
	28	C END
	29
	30	C-----------------------------------------------------------------------
	31
	32	C...SaSgam version 2 - parton distributions of the photon
	33
	34	C...by Gerhard A. Schuler and Torbjorn Sjostrand
	35
	36	C...For further information see Z. Phys. C68 (1995) 607
	37
	38	C...and CERN-TH/96-04 and LU TP 96-2.
	39
	40	C...Program last changed on 18 January 1996.
	41
	42	C
	43
	44	C!!!Note that one further call parameter - IP2 - has been added
	45
	46	C!!!to the SASGAM argument list compared with version 1.
	47
	48	C
	49
	50	C...The user should only need to call the SASGAM routine,
	51
	52	C...which in turn calls the auxiliary routines SASVMD, SASANO,
	53
	54	C...SASBEH and SASDIR. The package is self-contained.
	55
	56	C
	57
	58	C...One particular aspect of these parametrizations is that F2 for
	59
	60	C...the photon is not obtained just as the charge-squared-weighted
	61
	62	C...sum of quark distributions, but differ in the treatment of
	63
	64	C...heavy flavours (in F2 the DIS relation W2 = Q2*(1-x)/x restricts
65
66	C...the kinematics range of heavy-flavour production, but the same
67
68	C...kinematics is not relevant e.g. for jet production) and, for the
69
70	C...'MSbar' fits, in the addition of a Cgamma term related to the
71
72	C...separation of direct processes. Schematically:
73
74	C...PDF = VMD (rho, omega, phi) + anomalous (d, u, s, c, b).
75
76	C...F2 = VMD (rho, omega, phi) + anomalous (d, u, s) +
77
78	C... Bethe-Heitler (c, b) (+ Cgamma (d, u, s)).
79
80	C...The J/psi and Upsilon states have not been included in the VMD sum,
81
82	C...but low c and b masses in the other components should compensate
83
84	C...for this in a duality sense.
85
86	C
87
88	C...The calling sequence is the following:
89
90	C CALL SASGAM(ISET,X,Q2,P2,IP2,F2GM,XPDFGM)
91
92	C...with the following declaration statement:
93
94	C DIMENSION XPDFGM(-6:6)
95
96	C...and, optionally, further information in:
97
98	C COMMON/SASCOM/XPVMD(-6:6),XPANL(-6:6),XPANH(-6:6),XPBEH(-6:6),
99
100	C &XPDIR(-6:6)
101
102	C COMMON/SASVAL/VXPVMD(-6:6),VXPANL(-6:6),VXPANH(-6:6),VXPDGM(-6:6)
103
104	C...Input: ISET = 1 : SaS set 1D ('DIS', Q0 = 0.6 GeV)
105
106	C = 2 : SaS set 1M ('MSbar', Q0 = 0.6 GeV)
107
108	C = 3 : SaS set 2D ('DIS', Q0 = 2 GeV)
109
110	C = 4 : SaS set 2M ('MSbar', Q0 = 2 GeV)
111
112	C X : x value.
113
114	C Q2 : Q2 value.
115
116	C P2 : P2 value; should be = 0. for an on-shell photon.
117
118	C IP2 : scheme used to evaluate off-shell anomalous component.
119
120	C = 0 : recommended default, see = 7.
121
122	C = 1 : dipole dampening by integration; very time-consuming.
123
124	C = 2 : P_0^2 = max( Q_0^2, P^2 )
125
126	C = 3 : P'_0^2 = Q_0^2 + P^2.
127
128	C = 4 : P_{eff} that preserves momentum sum.
129
130	C = 5 : P_{int} that preserves momentum and average
131
132	C evolution range.
133
134	C = 6 : P_{eff}, matched to P_0 in P2 -> Q2 limit.
135
136	C = 7 : P_{eff}, matched to P_0 in P2 -> Q2 limit.
137
138	C...Output: F2GM : F2 value of the photon (including factors of alpha_em).
139
140	C XPFDGM : x times parton distribution functions of the photon,
141
142	C with elements 0 = g, 1 = d, 2 = u, 3 = s, 4 = c, 5 = b,
143
144	C 6 = t (always empty!), - for antiquarks (result is same).
145
146	C...The breakdown by component is stored in the commonblock SASCOM,
147
148	C with elements as above.
149
150	C XPVMD : rho, omega, phi VMD part only of output.
151
152	C XPANL : d, u, s anomalous part only of output.
153
154	C XPANH : c, b anomalous part only of output.
155
156	C XPBEH : c, b Bethe-Heitler part only of output.
157
158	C XPDIR : Cgamma (direct contribution) part only of output.
159
160	C...The above arrays do not distinguish valence and sea contributions,
161
162	C...although this information is available internally. The additional
163
164	C...commonblock SASVAL provides the valence part only of the above
165
166	C...distributions. Array names VXPVMD, VXPANL and VXPANH correspond
167
168	C...to XPVMD, XPANL and XPANH, while XPBEH and XPDIR are valence only
169
170	C...and therefore not given doubly. VXPDGM gives the sum of valence
171
172	C...parts, and so matches XPDFGM. The difference, i.e. XPVMD-VXPVMD
173
174	C...and so on, gives the sea part only.
175
176	C
177
178	SUBROUTINE SASGAM(ISET,X,Q2,P2,IP2,F2GM,XPDFGM)
179
180	C...Purpose: to construct the F2 and parton distributions of the photon
181
182	C...by summing homogeneous (VMD) and inhomogeneous (anomalous) terms.
183
184	C...For F2, c and b are included by the Bethe-Heitler formula;
185
186	C...in the 'MSbar' scheme additionally a Cgamma term is added.
187
188	DIMENSION XPDFGM(-6:6)
189
190	COMMON/SASCOM/XPVMD(-6:6),XPANL(-6:6),XPANH(-6:6),XPBEH(-6:6),
191
192	&XPDIR(-6:6)
193
194	COMMON/SASVAL/VXPVMD(-6:6),VXPANL(-6:6),VXPANH(-6:6),VXPDGM(-6:6)
195
196	SAVE /SASCOM/,/SASVAL/
197
198	C
199
200	C...Temporary array.
201
202	DIMENSION XPGA(-6:6), VXPGA(-6:6)
203
204	C...Charm and bottom masses (low to compensate for J/psi etc.).
205
206	DATA PMC/1.3/, PMB/4.6/
207
208	C...alpha_em and alpha_em/(2*pi).
209
210	DATA AEM/0.007297/, AEM2PI/0.0011614/
211
212	C...Lambda value for 4 flavours.
213
214	DATA ALAM/0.20/
215
216	C...Mixture u/(u+d), = 0.5 for incoherent and = 0.8 for coherent sum.
217
218	DATA FRACU/0.8/
219
220	C...VMD couplings f_V*2/(4pi).
221
222	DATA FRHO/2.20/, FOMEGA/23.6/, FPHI/18.4/
223
224	C...Masses for rho (=omega) and phi.
225
226	DATA PMRHO/0.770/, PMPHI/1.020/
227
228	C...Number of points in integration for IP2=1.
229
230	DATA NSTEP/100/
231
232	C
233
234	C...Reset output.
235
236	F2GM=0.
237
238	DO 100 KFL=-6,6
239
240	XPDFGM(KFL)=0.
241
242	XPVMD(KFL)=0.
243
244	XPANL(KFL)=0.
245
246	XPANH(KFL)=0.
247
248	XPBEH(KFL)=0.
249
250	XPDIR(KFL)=0.
251
252	VXPVMD(KFL)=0.
253
254	VXPANL(KFL)=0.
255
256	VXPANH(KFL)=0.
257
258	VXPDGM(KFL)=0.
259
260	100 CONTINUE
261
262	C
263
264	C...Check that input sensible.
265
266	IF(ISET.LE.0.OR.ISET.GE.5) THEN
267
268	WRITE(,) ' FATAL ERROR: SaSgam called for unknown set'
269
270	WRITE(,) ' ISET = ',ISET
271
272	STOP
273
274	ENDIF
275
276	IF(X.LE.0..OR.X.GT.1.) THEN
277
278	WRITE(,) ' FATAL ERROR: SaSgam called for unphysical x'
279
280	WRITE(,) ' X = ',X
281
282	STOP
283
284	ENDIF
285
286	C
287
288	C...Set Q0 cut-off parameter as function of set used.
289
290	IF(ISET.LE.2) THEN
291
292	Q0=0.6
293
294	ELSE
295
296	Q0=2.
297
298	ENDIF
299
300	Q02=Q0**2
301
302	C
303
304	C...Scale choice for off-shell photon; common factors.
305
306	Q2A=Q2
307
308	FACNOR=1.
309
310	IF(IP2.EQ.1) THEN
311
312	P2MX=P2+Q02
313
314	Q2A=Q2+P2*Q02/MAX(Q02,Q2)
315
316	FACNOR=LOG(Q2/Q02)/NSTEP
317
318	ELSEIF(IP2.EQ.2) THEN
319
320	P2MX=MAX(P2,Q02)
321
322	ELSEIF(IP2.EQ.3) THEN
323
324	P2MX=P2+Q02
325
326	Q2A=Q2+P2*Q02/MAX(Q02,Q2)
327
328	ELSEIF(IP2.EQ.4) THEN
329
330	P2MX=Q2(Q02+P2)/(Q2+P2)EXP(P2*(Q2-Q02)/
331
332	& ((Q2+P2)*(Q02+P2)))
333
334	ELSEIF(IP2.EQ.5) THEN
335
336	P2MXA=Q2(Q02+P2)/(Q2+P2)EXP(P2*(Q2-Q02)/
337
338	& ((Q2+P2)*(Q02+P2)))
339
340	P2MX=Q0*SQRT(P2MXA)
341
342	FACNOR=LOG(Q2/P2MXA)/LOG(Q2/P2MX)
343
344	ELSEIF(IP2.EQ.6) THEN
345
346	P2MX=Q2(Q02+P2)/(Q2+P2)EXP(P2*(Q2-Q02)/
347
348	& ((Q2+P2)*(Q02+P2)))
349
350	P2MX=MAX(0.,1.-P2/Q2)P2MX+MIN(1.,P2/Q2)MAX(P2,Q02)
351
352	ELSE
353
354	P2MXA=Q2(Q02+P2)/(Q2+P2)EXP(P2*(Q2-Q02)/
355
356	& ((Q2+P2)*(Q02+P2)))
357
358	P2MX=Q0*SQRT(P2MXA)
359
360	P2MXB=P2MX
361
362	P2MX=MAX(0.,1.-P2/Q2)P2MX+MIN(1.,P2/Q2)MAX(P2,Q02)
363
364	P2MXB=MAX(0.,1.-P2/Q2)P2MXB+MIN(1.,P2/Q2)P2MXA
365
366	FACNOR=LOG(Q2/P2MXA)/LOG(Q2/P2MXB)
367
368	ENDIF
369
370	C
371
372	C...Call VMD parametrization for d quark and use to give rho, omega,
373
374	C...phi. Note dipole dampening for off-shell photon.
375
376	CALL SASVMD(ISET,1,X,Q2A,P2MX,ALAM,XPGA,VXPGA)
377
378	XFVAL=VXPGA(1)
379
380	XPGA(1)=XPGA(2)
381
382	XPGA(-1)=XPGA(-2)
383
384	FACUD=AEM(1./FRHO+1./FOMEGA)(PMRHO2/(PMRHO2+P2))**2
385
386	FACS=AEM(1./FPHI)(PMPHI2/(PMPHI2+P2))**2
387
388	DO 110 KFL=-5,5
389
390	XPVMD(KFL)=(FACUD+FACS)*XPGA(KFL)
391
392	110 CONTINUE
393
394	XPVMD(1)=XPVMD(1)+(1.-FRACU)FACUDXFVAL
395
396	XPVMD(2)=XPVMD(2)+FRACUFACUDXFVAL
397
398	XPVMD(3)=XPVMD(3)+FACS*XFVAL
399
400	XPVMD(-1)=XPVMD(-1)+(1.-FRACU)FACUDXFVAL
401
402	XPVMD(-2)=XPVMD(-2)+FRACUFACUDXFVAL
403
404	XPVMD(-3)=XPVMD(-3)+FACS*XFVAL
405
406	VXPVMD(1)=(1.-FRACU)FACUDXFVAL
407
408	VXPVMD(2)=FRACUFACUDXFVAL
409
410	VXPVMD(3)=FACS*XFVAL
411
412	VXPVMD(-1)=(1.-FRACU)FACUDXFVAL
413
414	VXPVMD(-2)=FRACUFACUDXFVAL
415
416	VXPVMD(-3)=FACS*XFVAL
417
418	C
419
420	IF(IP2.NE.1) THEN
421
422	C...Anomalous parametrizations for different strategies
423
424	C...for off-shell photons; except full integration.
425
426	C
427
428	C...Call anomalous parametrization for d + u + s.
429
430	CALL SASANO(-3,X,Q2A,P2MX,ALAM,XPGA,VXPGA)
431
432	DO 120 KFL=-5,5
433
434	XPANL(KFL)=FACNOR*XPGA(KFL)
435
436	VXPANL(KFL)=FACNOR*VXPGA(KFL)
437
438	120 CONTINUE
439
440	C
441
442	C...Call anomalous parametrization for c and b.
443
444	CALL SASANO(4,X,Q2A,P2MX,ALAM,XPGA,VXPGA)
445
446	DO 130 KFL=-5,5
447
448	XPANH(KFL)=FACNOR*XPGA(KFL)
449
450	VXPANH(KFL)=FACNOR*VXPGA(KFL)
451
452	130 CONTINUE
453
454	CALL SASANO(5,X,Q2A,P2MX,ALAM,XPGA,VXPGA)
455
456	DO 140 KFL=-5,5
457
458	XPANH(KFL)=XPANH(KFL)+FACNOR*XPGA(KFL)
459
460	VXPANH(KFL)=VXPANH(KFL)+FACNOR*VXPGA(KFL)
461
462	140 CONTINUE
463
464	C
465
466	ELSE
467
468	C...Special option: loop over flavours and integrate over k2.
469
470	DO 170 KF=1,5
471
472	DO 160 ISTEP=1,NSTEP
473
474	Q2STEP=Q02(Q2/Q02)*((ISTEP-0.5)/NSTEP)
475
476	IF((KF.EQ.4.AND.Q2STEP.LT.PMC**2).OR.
477
478	& (KF.EQ.5.AND.Q2STEP.LT.PMB**2)) GOTO 160
479
480	CALL SASVMD(0,KF,X,Q2,Q2STEP,ALAM,XPGA,VXPGA)
481
482	FACQ=AEM2PI(Q2STEP/(Q2STEP+P2))2FACNOR
483
484	IF(MOD(KF,2).EQ.0) FACQ=FACQ*(8./9.)
485
486	IF(MOD(KF,2).EQ.1) FACQ=FACQ*(2./9.)
487
488	DO 150 KFL=-5,5
489
490	IF(KF.LE.3) XPANL(KFL)=XPANL(KFL)+FACQ*XPGA(KFL)
491
492	IF(KF.GE.4) XPANH(KFL)=XPANH(KFL)+FACQ*XPGA(KFL)
493
494	IF(KF.LE.3) VXPANL(KFL)=VXPANL(KFL)+FACQ*VXPGA(KFL)
495
496	IF(KF.GE.4) VXPANH(KFL)=VXPANH(KFL)+FACQ*VXPGA(KFL)
497
498	150 CONTINUE
499
500	160 CONTINUE
501
502	170 CONTINUE
503
504	ENDIF
505
506	C
507
508	C...Call Bethe-Heitler term expression for charm and bottom.
509
510	CALL SASBEH(4,X,Q2,P2,PMC**2,XPBH)
511
512	XPBEH(4)=XPBH
513
514	XPBEH(-4)=XPBH
515
516	CALL SASBEH(5,X,Q2,P2,PMB**2,XPBH)
517
518	XPBEH(5)=XPBH
519
520	XPBEH(-5)=XPBH
521
522	C
523
524	C...For MSbar subtraction call C^gamma term expression for d, u, s.
525
526	IF(ISET.EQ.2.OR.ISET.EQ.4) THEN
527
528	CALL SASDIR(X,Q2,P2,Q02,XPGA)
529
530	DO 180 KFL=-5,5
531
532	XPDIR(KFL)=XPGA(KFL)
533
534	180 CONTINUE
535
536	ENDIF
537
538	C
539
540	C...Store result in output array.
541
542	DO 190 KFL=-5,5
543
544	CHSQ=1./9.
545
546	IF(IABS(KFL).EQ.2.OR.IABS(KFL).EQ.4) CHSQ=4./9.
547
548	XPF2=XPVMD(KFL)+XPANL(KFL)+XPBEH(KFL)+XPDIR(KFL)
549
550	IF(KFL.NE.0) F2GM=F2GM+CHSQ*XPF2
551
552	XPDFGM(KFL)=XPVMD(KFL)+XPANL(KFL)+XPANH(KFL)
553
554	VXPDGM(KFL)=VXPVMD(KFL)+VXPANL(KFL)+VXPANH(KFL)
555
556	190 CONTINUE
557
558	C
559
560	RETURN
561
562	END