~maddevelopers/mg5amcnlo/WWW_tg

« back to all changes in this revision

Viewing changes to users/mardelcourt/PROC_242195/PROC_242195/Source/PDF/Ctq5Par.f

  • Committer: Martin Delcourt
  • Date: 2013-03-25 20:56:51 UTC
  • Revision ID: mdel@perso.be-20130325205651-f3bisj4z2jlfyjh3
added test et test2

Show diffs side-by-side

added added

removed removed

Lines of Context:
1
 
C   CTEQ5M1 and CTEQ5L Parton Distribution Functions in Parametrized Form
2
 
C                             
3
 
C               September 15, 1999
4
 
C
5
 
C   Ref: "GLOBAL QCD ANALYSIS OF PARTON STRUCTURE OF THE NUCLEON:
6
 
C         CTEQ5 PPARTON DISTRIBUTIONS"
7
 
C   hep-ph/9903282
8
 
C
9
 
C   The CTEQ5M1 set given here is an updated version of the original CTEQ5M
10
 
C     set posted, in the table version, on the Web page of CTEQ.
11
 
C     The differences between CTEQ5M and CTEQ5M1 are insignificant for almost
12
 
C     all applications. 
13
 
C   The improvement is in the QCD evolution which is now more accurate, and
14
 
C   which agrees completely with the benchmark work of the HERA 96/97 Workshop.
15
 
 
16
 
C   The differences between the parametrized and the corresponding table ver-
17
 
C sions (on which it is based) are of similar order as between the two version.
18
 
C    
19
 
C!! Because accurate parametrizations over a wide range of (x,Q) is hard to
20
 
C   obtain, only the most widely used sets CTEQ5M and CTEQ5L are available 
21
 
C   in parametrized form for now. 
22
 
 
23
 
C   These parametrizations were obtained by Jon Pumplin.
24
 
C
25
 
C                    ******************************
26
 
C  Iset   PDF        Description                 Alpha_s(Mz)  Lam4  Lam5
27
 
C ---------------------------------------------------------------------------
28
 
C   1    CTEQ5M1  Standard NLO MSbar scheme         0.118     326   226
29
 
C   3    CTEQ5L   Leading Order                     0.127     192   146
30
 
C ---------------------------------------------------------------------------
31
 
C   Note the Qcd-lambda values given for CTEQ5L is for the leading order
32
 
C     form of Alpha_s!!  Alpha_s(Mz) gives the absolute calibration.
33
 
 
34
 
C  The two Iset value are adopted to agree with the standard table versions.
35
 
 
36
 
C   The following user-callable routines are provided:
37
 
38
 
C     FUNCTION Ctq5Pd (Iset, Iprtn, X, Q, Irt) 
39
 
C         returns the PROBABILITY density for a GIVEN flavor;
40
 
C
41
 
C     FUNCTION Ctq5df (Iset, Iprtn, X, Q, Irt)
42
 
C         returns the MOMENTUM density of a GIVEN valence or sea distribution.
43
 
C
44
 
C     SUBROUTINE Ctq5Pds(Iset, Pdf, X, Q, Irt)
45
 
C         returns an array of MOMENTUM densities for ALL flavors;
46
 
C****************************************************************************
47
 
C     Added by Tim Stelzer    
48
 
c               subroutine pftopdg(x,Q,pdf)
49
 
c               subroutine pdfset(parm,val)
50
 
c
51
 
c
52
 
C
53
 
C   The arguments of these routines are as follows: 
54
 
C
55
 
C   Iset is the set number:  1 for CTEQ5M1 or 3 for CTEQ5L  
56
 
C
57
 
C   Iprtn  is the parton label (6, 5, 4, 3, 2, 1, 0, -1, ......, -6)
58
 
C                          for (t, b, c, s, d, u, g, u_bar, ..., t_bar)
59
 
C  *** WARNING: We use the parton label 2 as D-quark and 1 as U-quark, 
60
 
C               which might be different from your labels.
61
 
C
62
 
C   X, Q are the usual x, Q; 
63
 
C
64
 
C   Irt is an error code: 0 if there was no error; 1 or more if (x,q) was 
65
 
C   outside the range of validity of the parametrization.
66
 
C       
67
 
C  Range of validity:
68
 
C  
69
 
C     The range of (x, Q) covered by this parametrization of the QCD evolved
70
 
C     parton distributions is 1E-6 < x < 1 ; 1.1 GeV < Q < 10 TeV.  Of course,
71
 
C     the PDF's are constrained by data only in a subset of that region; and 
72
 
C     the assumed DGLAP evolution is unlikely to be valid for all of it either.
73
 
C
74
 
C     The range of (x, Q) used in the CTEQ5 round of global analysis is 
75
 
C     approximately 0.01 < x < 0.75 ; and 4 GeV^2 < Q^2 < 400 GeV^2 for 
76
 
C     fixed target experiments; 0.0001 < x < 0.3 from HERA data; and   
77
 
C     Q^2 up to 40,000 GeV^2 from Tevatron inclusive Jet data.
78
 
C
79
 
C   DOUBLE PRECISION is used throughout in these routines, but conversion to 
80
 
C   SINGLE PRECISION is possible by removing the Implicit Double Precision statements. 
81
 
C
82
 
C **************************************************************************
83
 
 
84
 
 
85
 
C ********************************************************
86
 
      FUNCTION CTQ5PD(ISET, IPARTON, X, Q, IRT)
87
 
C ********************************************************
88
 
      IMPLICIT DOUBLE PRECISION (A-H, O-Z)
89
 
 
90
 
c if called at a point (x,q) that is outside the region that was 
91
 
c actually parametrized, return a value of 0, and set the error code IRT=1.  
92
 
c The user can remove the following IF statement to receive instead an 
93
 
c extrapolated value, which may be wildly unphysical.
94
 
      if((x .lt. 1.e-6). or. (x .gt. 1.) 
95
 
     &   .or. (q .lt. .99) .or. (q .gt. 10000.)) then
96
 
         ctq5pd = 0.d0
97
 
         irt = 1
98
 
         return
99
 
      endif
100
 
 
101
 
      irt = 0
102
 
      if(iset .eq. 3) then
103
 
         ctq5pd = ctq5L(iparton,x,q)
104
 
      elseif(iset .eq. 1) then
105
 
         ctq5pd = ctq5Mi(iparton,x,q)
106
 
      else
107
 
         print *,'iset=',iset,' has not been parametrized.' 
108
 
           print '(/A)', 'Use the interpolation-table version instead.'
109
 
         stop
110
 
      endif
111
 
 
112
 
      return
113
 
      end
114
 
 
115
 
C ********************************************************
116
 
      FUNCTION CTQ5DF(ISET, IFL, X, Q, IRT)
117
 
C ********************************************************
118
 
      IMPLICIT DOUBLE PRECISION (A-H, O-Z)
119
 
 
120
 
      CTQ5DF = X * CTQ5PD(ISET, IPARTON, X, Q, IRT)
121
 
        
122
 
      RETURN
123
 
      END
124
 
 
125
 
C ********************************************************
126
 
      SUBROUTINE CTQ5PDS(ISET, PDF, X, Q, IRT)
127
 
C ********************************************************
128
 
      IMPLICIT DOUBLE PRECISION (A-H, O-Z)
129
 
      DIMENSION PDF (-6:6)
130
 
 
131
 
      IRT = 0
132
 
 
133
 
      DO IFL= -6,2
134
 
         PDF(IFL) = CTQ5PD(ISET,IFL,X,Q,IRT1)
135
 
         IRT = IRT + IRT1
136
 
 
137
 
         IF (IFL .LE. -3) THEN
138
 
            PDF(-IFL) = PDF(IFL)
139
 
         ENDIF
140
 
 
141
 
      ENDDO
142
 
 
143
 
      RETURN
144
 
      END
145
 
 
146
 
c --------------------------------------------------------------------------
147
 
        double precision function ctq5MI(ifl,x,q)
148
 
c Parametrization of cteq5MI parton distribution functions (J. Pumplin 9/99).
149
 
c ifl: 1=u,2=d,3=s,4=c,5=b;0=gluon;-1=ubar,-2=dbar,-3=sbar,-4=cbar,-5=bbar.
150
 
c --------------------------------------------------------------------------
151
 
        implicit double precision (a-h,o-z)
152
 
        integer ifl
153
 
 
154
 
        ii = ifl
155
 
        if(ii .gt. 2) then
156
 
           ii = -ii
157
 
        endif
158
 
 
159
 
        if(ii .eq. -1) then
160
 
           sum = faux5MI(-1,x,q)
161
 
           ratio = faux5MI(-2,x,q)
162
 
           ctq5MI = sum/(1.d0 + ratio)
163
 
 
164
 
        elseif(ii .eq. -2) then
165
 
           sum = faux5MI(-1,x,q)
166
 
           ratio = faux5MI(-2,x,q)
167
 
           ctq5MI = sum*ratio/(1.d0 + ratio)
168
 
 
169
 
        elseif(ii .ge. -5) then
170
 
           ctq5MI = faux5MI(ii,x,q)
171
 
 
172
 
        else
173
 
           ctq5MI = 0.d0 
174
 
 
175
 
        endif
176
 
 
177
 
        return
178
 
        end
179
 
 
180
 
c ---------------------------------------------------------------------
181
 
      double precision function faux5MI(ifl,x,q)
182
 
c auxiliary function for parametrization of CTEQ5MI (J. Pumplin 9/99).
183
 
c ---------------------------------------------------------------------
184
 
      implicit double precision (a-h,o-z)
185
 
      integer ifl
186
 
 
187
 
      parameter (nex=8, nlf=2)
188
 
      dimension am(0:nex,0:nlf,-5:2)
189
 
      dimension alfvec(-5:2), qmavec(-5:2)
190
 
      dimension mexvec(-5:2), mlfvec(-5:2)
191
 
      dimension ut1vec(-5:2), ut2vec(-5:2)
192
 
      dimension af(0:nex)
193
 
 
194
 
      data mexvec( 2) / 8 /
195
 
      data mlfvec( 2) / 2 /
196
 
      data ut1vec( 2) /  0.5141718E+01 /
197
 
      data ut2vec( 2) / -0.1346944E+01 /
198
 
      data alfvec( 2) /  0.5260555E+00 /
199
 
      data qmavec( 2) /  0.0000000E+00 /
200
 
      data (am( 0,k, 2),k=0, 2)
201
 
     & /  0.4289071E+01, -0.2536870E+01, -0.1259948E+01 /
202
 
      data (am( 1,k, 2),k=0, 2)
203
 
     & /  0.9839410E+00,  0.4168426E-01, -0.5018952E-01 /
204
 
      data (am( 2,k, 2),k=0, 2)
205
 
     & / -0.1651961E+02,  0.9246261E+01,  0.5996400E+01 /
206
 
      data (am( 3,k, 2),k=0, 2)
207
 
     & / -0.2077936E+02,  0.9786469E+01,  0.7656465E+01 /
208
 
      data (am( 4,k, 2),k=0, 2)
209
 
     & /  0.3054926E+02,  0.1889536E+01,  0.1380541E+01 /
210
 
      data (am( 5,k, 2),k=0, 2)
211
 
     & /  0.3084695E+02, -0.1212303E+02, -0.1053551E+02 /
212
 
      data (am( 6,k, 2),k=0, 2)
213
 
     & / -0.1426778E+02,  0.6239537E+01,  0.5254819E+01 /
214
 
      data (am( 7,k, 2),k=0, 2)
215
 
     & / -0.1909811E+02,  0.3695678E+01,  0.5495729E+01 /
216
 
      data (am( 8,k, 2),k=0, 2)
217
 
     & /  0.1889751E-01,  0.5027193E-02,  0.6624896E-03 /
218
 
 
219
 
      data mexvec( 1) / 8 /
220
 
      data mlfvec( 1) / 2 /
221
 
      data ut1vec( 1) /  0.4138426E+01 /
222
 
      data ut2vec( 1) / -0.3221374E+01 /
223
 
      data alfvec( 1) /  0.4960962E+00 /
224
 
      data qmavec( 1) /  0.0000000E+00 /
225
 
      data (am( 0,k, 1),k=0, 2)
226
 
     & /  0.1332497E+01, -0.3703718E+00,  0.1288638E+00 /
227
 
      data (am( 1,k, 1),k=0, 2)
228
 
     & /  0.7544687E+00,  0.3255075E-01, -0.4706680E-01 /
229
 
      data (am( 2,k, 1),k=0, 2)
230
 
     & / -0.7638814E+00,  0.5008313E+00, -0.9237374E-01 /
231
 
      data (am( 3,k, 1),k=0, 2)
232
 
     & / -0.3689889E+00, -0.1055098E+01, -0.4645065E+00 /
233
 
      data (am( 4,k, 1),k=0, 2)
234
 
     & /  0.3991610E+02,  0.1979881E+01,  0.1775814E+01 /
235
 
      data (am( 5,k, 1),k=0, 2)
236
 
     & /  0.6201080E+01,  0.2046288E+01,  0.3804571E+00 /
237
 
      data (am( 6,k, 1),k=0, 2)
238
 
     & / -0.8027900E+00, -0.7011688E+00, -0.8049612E+00 /
239
 
      data (am( 7,k, 1),k=0, 2)
240
 
     & / -0.8631305E+01, -0.3981200E+01,  0.6970153E+00 /
241
 
      data (am( 8,k, 1),k=0, 2)
242
 
     & /  0.2371230E-01,  0.5372683E-02,  0.1118701E-02 /
243
 
 
244
 
      data mexvec( 0) / 8 /
245
 
      data mlfvec( 0) / 2 /
246
 
      data ut1vec( 0) / -0.1026789E+01 /
247
 
      data ut2vec( 0) / -0.9051707E+01 /
248
 
      data alfvec( 0) /  0.9462977E+00 /
249
 
      data qmavec( 0) /  0.0000000E+00 /
250
 
      data (am( 0,k, 0),k=0, 2)
251
 
     & /  0.1191990E+03, -0.8548739E+00, -0.1963040E+01 /
252
 
      data (am( 1,k, 0),k=0, 2)
253
 
     & / -0.9449972E+02,  0.1074771E+01,  0.2056055E+01 /
254
 
      data (am( 2,k, 0),k=0, 2)
255
 
     & /  0.3701064E+01, -0.1167947E-02,  0.1933573E+00 /
256
 
      data (am( 3,k, 0),k=0, 2)
257
 
     & /  0.1171345E+03, -0.1064540E+01, -0.1875312E+01 /
258
 
      data (am( 4,k, 0),k=0, 2)
259
 
     & / -0.1014453E+03, -0.5707427E+00,  0.4511242E-01 /
260
 
      data (am( 5,k, 0),k=0, 2)
261
 
     & /  0.6365168E+01,  0.1275354E+01, -0.4964081E+00 /
262
 
      data (am( 6,k, 0),k=0, 2)
263
 
     & / -0.3370693E+01, -0.1122020E+01,  0.5947751E-01 /
264
 
      data (am( 7,k, 0),k=0, 2)
265
 
     & / -0.5327270E+01, -0.9293556E+00,  0.6629940E+00 /
266
 
      data (am( 8,k, 0),k=0, 2)
267
 
     & /  0.2437513E-01,  0.1600939E-02,  0.6855336E-03 /
268
 
 
269
 
      data mexvec(-1) / 8 /
270
 
      data mlfvec(-1) / 2 /
271
 
      data ut1vec(-1) /  0.5243571E+01 /
272
 
      data ut2vec(-1) / -0.2870513E+01 /
273
 
      data alfvec(-1) /  0.6701448E+00 /
274
 
      data qmavec(-1) /  0.0000000E+00 /
275
 
      data (am( 0,k,-1),k=0, 2)
276
 
     & /  0.2428863E+02,  0.1907035E+01, -0.4606457E+00 /
277
 
      data (am( 1,k,-1),k=0, 2)
278
 
     & /  0.2006810E+01, -0.1265915E+00,  0.7153556E-02 /
279
 
      data (am( 2,k,-1),k=0, 2)
280
 
     & / -0.1884546E+02, -0.2339471E+01,  0.5740679E+01 /
281
 
      data (am( 3,k,-1),k=0, 2)
282
 
     & / -0.2527892E+02, -0.2044124E+01,  0.1280470E+02 /
283
 
      data (am( 4,k,-1),k=0, 2)
284
 
     & / -0.1013824E+03, -0.1594199E+01,  0.2216401E+00 /
285
 
      data (am( 5,k,-1),k=0, 2)
286
 
     & /  0.8070930E+02,  0.1792072E+01, -0.2164364E+02 /
287
 
      data (am( 6,k,-1),k=0, 2)
288
 
     & / -0.4641050E+02,  0.1977338E+00,  0.1273014E+02 /
289
 
      data (am( 7,k,-1),k=0, 2)
290
 
     & / -0.3910568E+02,  0.1719632E+01,  0.1086525E+02 /
291
 
      data (am( 8,k,-1),k=0, 2)
292
 
     & / -0.1185496E+01, -0.1905847E+00, -0.8744118E-03 /
293
 
 
294
 
      data mexvec(-2) / 7 /
295
 
      data mlfvec(-2) / 2 /
296
 
      data ut1vec(-2) /  0.4782210E+01 /
297
 
      data ut2vec(-2) / -0.1976856E+02 /
298
 
      data alfvec(-2) /  0.7558374E+00 /
299
 
      data qmavec(-2) /  0.0000000E+00 /
300
 
      data (am( 0,k,-2),k=0, 2)
301
 
     & / -0.6216935E+00,  0.2369963E+00, -0.7909949E-02 /
302
 
      data (am( 1,k,-2),k=0, 2)
303
 
     & /  0.1245440E+01, -0.1031510E+00,  0.4916523E-02 /
304
 
      data (am( 2,k,-2),k=0, 2)
305
 
     & / -0.7060824E+01, -0.3875283E-01,  0.1784981E+00 /
306
 
      data (am( 3,k,-2),k=0, 2)
307
 
     & / -0.7430595E+01,  0.1964572E+00, -0.1284999E+00 /
308
 
      data (am( 4,k,-2),k=0, 2)
309
 
     & / -0.6897810E+01,  0.2620543E+01,  0.8012553E-02 /
310
 
      data (am( 5,k,-2),k=0, 2)
311
 
     & /  0.1507713E+02,  0.2340307E-01,  0.2482535E+01 /
312
 
      data (am( 6,k,-2),k=0, 2)
313
 
     & / -0.1815341E+01, -0.1538698E+01, -0.2014208E+01 /
314
 
      data (am( 7,k,-2),k=0, 2)
315
 
     & / -0.2571932E+02,  0.2903941E+00, -0.2848206E+01 /
316
 
 
317
 
      data mexvec(-3) / 7 /
318
 
      data mlfvec(-3) / 2 /
319
 
      data ut1vec(-3) /  0.4518239E+01 /
320
 
      data ut2vec(-3) / -0.2690590E+01 /
321
 
      data alfvec(-3) /  0.6124079E+00 /
322
 
      data qmavec(-3) /  0.0000000E+00 /
323
 
      data (am( 0,k,-3),k=0, 2)
324
 
     & / -0.2734458E+01, -0.7245673E+00, -0.6351374E+00 /
325
 
      data (am( 1,k,-3),k=0, 2)
326
 
     & /  0.2927174E+01,  0.4822709E+00, -0.1088787E-01 /
327
 
      data (am( 2,k,-3),k=0, 2)
328
 
     & / -0.1771017E+02, -0.1416635E+01,  0.8467622E+01 /
329
 
      data (am( 3,k,-3),k=0, 2)
330
 
     & / -0.4972782E+02, -0.3348547E+01,  0.1767061E+02 /
331
 
      data (am( 4,k,-3),k=0, 2)
332
 
     & / -0.7102770E+01, -0.3205337E+01,  0.4101704E+00 /
333
 
      data (am( 5,k,-3),k=0, 2)
334
 
     & /  0.7169698E+02, -0.2205985E+01, -0.2463931E+02 /
335
 
      data (am( 6,k,-3),k=0, 2)
336
 
     & / -0.4090347E+02,  0.2103486E+01,  0.1416507E+02 /
337
 
      data (am( 7,k,-3),k=0, 2)
338
 
     & / -0.2952639E+02,  0.5376136E+01,  0.7825585E+01 /
339
 
 
340
 
      data mexvec(-4) / 7 /
341
 
      data mlfvec(-4) / 2 /
342
 
      data ut1vec(-4) /  0.2783230E+01 /
343
 
      data ut2vec(-4) / -0.1746328E+01 /
344
 
      data alfvec(-4) /  0.1115653E+01 /
345
 
      data qmavec(-4) /  0.1300000E+01 /
346
 
      data (am( 0,k,-4),k=0, 2)
347
 
     & / -0.1743872E+01, -0.1128921E+01, -0.2841969E+00 /
348
 
      data (am( 1,k,-4),k=0, 2)
349
 
     & /  0.3345755E+01,  0.3187765E+00,  0.1378124E+00 /
350
 
      data (am( 2,k,-4),k=0, 2)
351
 
     & / -0.2037615E+02,  0.4121687E+01,  0.2236520E+00 /
352
 
      data (am( 3,k,-4),k=0, 2)
353
 
     & / -0.4703104E+02,  0.5353087E+01, -0.1455347E+01 /
354
 
      data (am( 4,k,-4),k=0, 2)
355
 
     & / -0.1060230E+02, -0.1551122E+01, -0.1078863E+01 /
356
 
      data (am( 5,k,-4),k=0, 2)
357
 
     & /  0.5088892E+02, -0.8197304E+01,  0.8083451E+01 /
358
 
      data (am( 6,k,-4),k=0, 2)
359
 
     & / -0.2819070E+02,  0.4554086E+01, -0.5890995E+01 /
360
 
      data (am( 7,k,-4),k=0, 2)
361
 
     & / -0.1098238E+02,  0.2590096E+01, -0.8062879E+01 /
362
 
 
363
 
      data mexvec(-5) / 6 /
364
 
      data mlfvec(-5) / 2 /
365
 
      data ut1vec(-5) /  0.1619654E+02 /
366
 
      data ut2vec(-5) / -0.3367346E+01 /
367
 
      data alfvec(-5) /  0.5109891E-02 /
368
 
      data qmavec(-5) /  0.4500000E+01 /
369
 
      data (am( 0,k,-5),k=0, 2)
370
 
     & / -0.6800138E+01,  0.2493627E+01, -0.1075724E+01 /
371
 
      data (am( 1,k,-5),k=0, 2)
372
 
     & /  0.3036555E+01,  0.3324733E+00,  0.2008298E+00 /
373
 
      data (am( 2,k,-5),k=0, 2)
374
 
     & / -0.5203879E+01, -0.8493476E+01, -0.4523208E+01 /
375
 
      data (am( 3,k,-5),k=0, 2)
376
 
     & / -0.1524239E+01, -0.3411912E+01, -0.1771867E+02 /
377
 
      data (am( 4,k,-5),k=0, 2)
378
 
     & / -0.1099444E+02,  0.1320930E+01, -0.2353831E+01 /
379
 
      data (am( 5,k,-5),k=0, 2)
380
 
     & /  0.1699299E+02, -0.3565802E+02,  0.3566872E+02 /
381
 
      data (am( 6,k,-5),k=0, 2)
382
 
     & / -0.1465793E+02,  0.2703365E+02, -0.2176372E+02 /
383
 
 
384
 
      if(q .le. qmavec(ifl)) then
385
 
         faux5MI = 0.d0
386
 
         return
387
 
      endif
388
 
 
389
 
      if(x .ge. 1.d0) then
390
 
         faux5MI = 0.d0
391
 
         return
392
 
      endif
393
 
 
394
 
      tmp = log(q/alfvec(ifl))
395
 
      if(tmp .le. 0.d0) then
396
 
         faux5MI = 0.d0
397
 
         return
398
 
      endif
399
 
 
400
 
      sb = log(tmp)
401
 
      sb1 = sb - 1.2d0
402
 
      sb2 = sb1*sb1
403
 
 
404
 
      do i = 0, nex
405
 
         af(i) = 0.d0
406
 
         sbx = 1.d0
407
 
         do k = 0, mlfvec(ifl)
408
 
            af(i) = af(i) + sbx*am(i,k,ifl)
409
 
            sbx = sb1*sbx
410
 
         enddo
411
 
      enddo
412
 
 
413
 
      y = -log(x)
414
 
      u = log(x/0.00001d0)
415
 
 
416
 
      part1 = af(1)*y**(1.d0+0.01d0*af(4))*(1.d0+ af(8)*u)
417
 
      part2 = af(0)*(1.d0 - x) + af(3)*x 
418
 
      part3 = x*(1.d0-x)*(af(5)+af(6)*(1.d0-x)+af(7)*x*(1.d0-x))
419
 
      part4 = ut1vec(ifl)*log(1.d0-x) + 
420
 
     &        AF(2)*log(1.d0+exp(ut2vec(ifl))-x)
421
 
 
422
 
      faux5MI = exp(log(x) + part1 + part2 + part3 + part4)
423
 
 
424
 
c include threshold factor...
425
 
      faux5MI = faux5MI * (1.d0 - qmavec(ifl)/q)
426
 
 
427
 
      return
428
 
      end
429
 
c --------------------------------------------------------------------------
430
 
        double precision function ctq5L(ifl,x,q)
431
 
c Parametrization of cteq5L parton distribution functions (J. Pumplin 9/99).
432
 
c ifl: 1=u,2=d,3=s,4=c,5=b;0=gluon;-1=ubar,-2=dbar,-3=sbar,-4=cbar,-5=bbar.
433
 
c --------------------------------------------------------------------------
434
 
        implicit double precision (a-h,o-z)
435
 
        integer ifl
436
 
 
437
 
        ii = ifl
438
 
        if(ii .gt. 2) then
439
 
           ii = -ii
440
 
        endif
441
 
 
442
 
        if(ii .eq. -1) then
443
 
           sum = faux5L(-1,x,q)
444
 
           ratio = faux5L(-2,x,q)
445
 
           ctq5L = sum/(1.d0 + ratio)
446
 
 
447
 
        elseif(ii .eq. -2) then
448
 
           sum = faux5L(-1,x,q)
449
 
           ratio = faux5L(-2,x,q)
450
 
           ctq5L = sum*ratio/(1.d0 + ratio)
451
 
 
452
 
        elseif(ii .ge. -5) then
453
 
           ctq5L = faux5L(ii,x,q)
454
 
 
455
 
        else
456
 
           ctq5L = 0.d0 
457
 
 
458
 
        endif
459
 
 
460
 
        return
461
 
        end
462
 
 
463
 
c ---------------------------------------------------------------------
464
 
      double precision function faux5L(ifl,x,q)
465
 
c auxiliary function for parametrization of CTEQ5L (J. Pumplin 9/99).
466
 
c ---------------------------------------------------------------------
467
 
      implicit double precision (a-h,o-z)
468
 
      integer ifl
469
 
 
470
 
      parameter (nex=8, nlf=2)
471
 
      dimension am(0:nex,0:nlf,-5:2)
472
 
      dimension alfvec(-5:2), qmavec(-5:2)
473
 
      dimension mexvec(-5:2), mlfvec(-5:2)
474
 
      dimension ut1vec(-5:2), ut2vec(-5:2)
475
 
      dimension af(0:nex)
476
 
 
477
 
      data mexvec( 2) / 8 /
478
 
      data mlfvec( 2) / 2 /
479
 
      data ut1vec( 2) /  0.4971265E+01 /
480
 
      data ut2vec( 2) / -0.1105128E+01 /
481
 
      data alfvec( 2) /  0.2987216E+00 /
482
 
      data qmavec( 2) /  0.0000000E+00 /
483
 
      data (am( 0,k, 2),k=0, 2)
484
 
     & /  0.5292616E+01, -0.2751910E+01, -0.2488990E+01 /
485
 
      data (am( 1,k, 2),k=0, 2)
486
 
     & /  0.9714424E+00,  0.1011827E-01, -0.1023660E-01 /
487
 
      data (am( 2,k, 2),k=0, 2)
488
 
     & / -0.1651006E+02,  0.7959721E+01,  0.8810563E+01 /
489
 
      data (am( 3,k, 2),k=0, 2)
490
 
     & / -0.1643394E+02,  0.5892854E+01,  0.9348874E+01 /
491
 
      data (am( 4,k, 2),k=0, 2)
492
 
     & /  0.3067422E+02,  0.4235796E+01, -0.5112136E+00 /
493
 
      data (am( 5,k, 2),k=0, 2)
494
 
     & /  0.2352526E+02, -0.5305168E+01, -0.1169174E+02 /
495
 
      data (am( 6,k, 2),k=0, 2)
496
 
     & / -0.1095451E+02,  0.3006577E+01,  0.5638136E+01 /
497
 
      data (am( 7,k, 2),k=0, 2)
498
 
     & / -0.1172251E+02, -0.2183624E+01,  0.4955794E+01 /
499
 
      data (am( 8,k, 2),k=0, 2)
500
 
     & /  0.1662533E-01,  0.7622870E-02, -0.4895887E-03 /
501
 
 
502
 
      data mexvec( 1) / 8 /
503
 
      data mlfvec( 1) / 2 /
504
 
      data ut1vec( 1) /  0.2612618E+01 /
505
 
      data ut2vec( 1) / -0.1258304E+06 /
506
 
      data alfvec( 1) /  0.3407552E+00 /
507
 
      data qmavec( 1) /  0.0000000E+00 /
508
 
      data (am( 0,k, 1),k=0, 2)
509
 
     & /  0.9905300E+00, -0.4502235E+00,  0.1624441E+00 /
510
 
      data (am( 1,k, 1),k=0, 2)
511
 
     & /  0.8867534E+00,  0.1630829E-01, -0.4049085E-01 /
512
 
      data (am( 2,k, 1),k=0, 2)
513
 
     & /  0.8547974E+00,  0.3336301E+00,  0.1371388E+00 /
514
 
      data (am( 3,k, 1),k=0, 2)
515
 
     & /  0.2941113E+00, -0.1527905E+01,  0.2331879E+00 /
516
 
      data (am( 4,k, 1),k=0, 2)
517
 
     & /  0.3384235E+02,  0.3715315E+01,  0.8276930E+00 /
518
 
      data (am( 5,k, 1),k=0, 2)
519
 
     & /  0.6230115E+01,  0.3134639E+01, -0.1729099E+01 /
520
 
      data (am( 6,k, 1),k=0, 2)
521
 
     & / -0.1186928E+01, -0.3282460E+00,  0.1052020E+00 /
522
 
      data (am( 7,k, 1),k=0, 2)
523
 
     & / -0.8545702E+01, -0.6247947E+01,  0.3692561E+01 /
524
 
      data (am( 8,k, 1),k=0, 2)
525
 
     & /  0.1724598E-01,  0.7120465E-02,  0.4003646E-04 /
526
 
 
527
 
      data mexvec( 0) / 8 /
528
 
      data mlfvec( 0) / 2 /
529
 
      data ut1vec( 0) / -0.4656819E+00 /
530
 
      data ut2vec( 0) / -0.2742390E+03 /
531
 
      data alfvec( 0) /  0.4491863E+00 /
532
 
      data qmavec( 0) /  0.0000000E+00 /
533
 
      data (am( 0,k, 0),k=0, 2)
534
 
     & /  0.1193572E+03, -0.3886845E+01, -0.1133965E+01 /
535
 
      data (am( 1,k, 0),k=0, 2)
536
 
     & / -0.9421449E+02,  0.3995885E+01,  0.1607363E+01 /
537
 
      data (am( 2,k, 0),k=0, 2)
538
 
     & /  0.4206383E+01,  0.2485954E+00,  0.2497468E+00 /
539
 
      data (am( 3,k, 0),k=0, 2)
540
 
     & /  0.1210557E+03, -0.3015765E+01, -0.1423651E+01 /
541
 
      data (am( 4,k, 0),k=0, 2)
542
 
     & / -0.1013897E+03, -0.7113478E+00,  0.2621865E+00 /
543
 
      data (am( 5,k, 0),k=0, 2)
544
 
     & / -0.1312404E+01, -0.9297691E+00, -0.1562531E+00 /
545
 
      data (am( 6,k, 0),k=0, 2)
546
 
     & /  0.1627137E+01,  0.4954111E+00, -0.6387009E+00 /
547
 
      data (am( 7,k, 0),k=0, 2)
548
 
     & /  0.1537698E+00, -0.2487878E+00,  0.8305947E+00 /
549
 
      data (am( 8,k, 0),k=0, 2)
550
 
     & /  0.2496448E-01,  0.2457823E-02,  0.8234276E-03 /
551
 
 
552
 
      data mexvec(-1) / 8 /
553
 
      data mlfvec(-1) / 2 /
554
 
      data ut1vec(-1) /  0.3862583E+01 /
555
 
      data ut2vec(-1) / -0.1265969E+01 /
556
 
      data alfvec(-1) /  0.2457668E+00 /
557
 
      data qmavec(-1) /  0.0000000E+00 /
558
 
      data (am( 0,k,-1),k=0, 2)
559
 
     & /  0.2647441E+02,  0.1059277E+02, -0.9176654E+00 /
560
 
      data (am( 1,k,-1),k=0, 2)
561
 
     & /  0.1990636E+01,  0.8558918E-01,  0.4248667E-01 /
562
 
      data (am( 2,k,-1),k=0, 2)
563
 
     & / -0.1476095E+02, -0.3276255E+02,  0.1558110E+01 /
564
 
      data (am( 3,k,-1),k=0, 2)
565
 
     & / -0.2966889E+01, -0.3649037E+02,  0.1195914E+01 /
566
 
      data (am( 4,k,-1),k=0, 2)
567
 
     & / -0.1000519E+03, -0.2464635E+01,  0.1964849E+00 /
568
 
      data (am( 5,k,-1),k=0, 2)
569
 
     & /  0.3718331E+02,  0.4700389E+02, -0.2772142E+01 /
570
 
      data (am( 6,k,-1),k=0, 2)
571
 
     & / -0.1872722E+02, -0.2291189E+02,  0.1089052E+01 /
572
 
      data (am( 7,k,-1),k=0, 2)
573
 
     & / -0.1628146E+02, -0.1823993E+02,  0.2537369E+01 /
574
 
      data (am( 8,k,-1),k=0, 2)
575
 
     & / -0.1156300E+01, -0.1280495E+00,  0.5153245E-01 /
576
 
 
577
 
      data mexvec(-2) / 7 /
578
 
      data mlfvec(-2) / 2 /
579
 
      data ut1vec(-2) /  0.1895615E+00 /
580
 
      data ut2vec(-2) / -0.3069097E+01 /
581
 
      data alfvec(-2) /  0.5293999E+00 /
582
 
      data qmavec(-2) /  0.0000000E+00 /
583
 
      data (am( 0,k,-2),k=0, 2)
584
 
     & / -0.6556775E+00,  0.2490190E+00,  0.3966485E-01 /
585
 
      data (am( 1,k,-2),k=0, 2)
586
 
     & /  0.1305102E+01, -0.1188925E+00, -0.4600870E-02 /
587
 
      data (am( 2,k,-2),k=0, 2)
588
 
     & / -0.2371436E+01,  0.3566814E+00, -0.2834683E+00 /
589
 
      data (am( 3,k,-2),k=0, 2)
590
 
     & / -0.6152826E+01,  0.8339877E+00, -0.7233230E+00 /
591
 
      data (am( 4,k,-2),k=0, 2)
592
 
     & / -0.8346558E+01,  0.2892168E+01,  0.2137099E+00 /
593
 
      data (am( 5,k,-2),k=0, 2)
594
 
     & /  0.1279530E+02,  0.1021114E+00,  0.5787439E+00 /
595
 
      data (am( 6,k,-2),k=0, 2)
596
 
     & /  0.5858816E+00, -0.1940375E+01, -0.4029269E+00 /
597
 
      data (am( 7,k,-2),k=0, 2)
598
 
     & / -0.2795725E+02, -0.5263392E+00,  0.1290229E+01 /
599
 
 
600
 
      data mexvec(-3) / 7 /
601
 
      data mlfvec(-3) / 2 /
602
 
      data ut1vec(-3) /  0.3753257E+01 /
603
 
      data ut2vec(-3) / -0.1113085E+01 /
604
 
      data alfvec(-3) /  0.3713141E+00 /
605
 
      data qmavec(-3) /  0.0000000E+00 /
606
 
      data (am( 0,k,-3),k=0, 2)
607
 
     & /  0.1580931E+01, -0.2273826E+01, -0.1822245E+01 /
608
 
      data (am( 1,k,-3),k=0, 2)
609
 
     & /  0.2702644E+01,  0.6763243E+00,  0.7231586E-02 /
610
 
      data (am( 2,k,-3),k=0, 2)
611
 
     & / -0.1857924E+02,  0.3907500E+01,  0.5850109E+01 /
612
 
      data (am( 3,k,-3),k=0, 2)
613
 
     & / -0.3044793E+02,  0.2639332E+01,  0.5566644E+01 /
614
 
      data (am( 4,k,-3),k=0, 2)
615
 
     & / -0.4258011E+01, -0.5429244E+01,  0.4418946E+00 /
616
 
      data (am( 5,k,-3),k=0, 2)
617
 
     & /  0.3465259E+02, -0.5532604E+01, -0.4904153E+01 /
618
 
      data (am( 6,k,-3),k=0, 2)
619
 
     & / -0.1658858E+02,  0.2923275E+01,  0.2266286E+01 /
620
 
      data (am( 7,k,-3),k=0, 2)
621
 
     & / -0.1149263E+02,  0.2877475E+01, -0.7999105E+00 /
622
 
 
623
 
      data mexvec(-4) / 7 /
624
 
      data mlfvec(-4) / 2 /
625
 
      data ut1vec(-4) /  0.4400772E+01 /
626
 
      data ut2vec(-4) / -0.1356116E+01 /
627
 
      data alfvec(-4) /  0.3712017E-01 /
628
 
      data qmavec(-4) /  0.1300000E+01 /
629
 
      data (am( 0,k,-4),k=0, 2)
630
 
     & / -0.8293661E+00, -0.3982375E+01, -0.6494283E-01 /
631
 
      data (am( 1,k,-4),k=0, 2)
632
 
     & /  0.2754618E+01,  0.8338636E+00, -0.6885160E-01 /
633
 
      data (am( 2,k,-4),k=0, 2)
634
 
     & / -0.1657987E+02,  0.1439143E+02, -0.6887240E+00 /
635
 
      data (am( 3,k,-4),k=0, 2)
636
 
     & / -0.2800703E+02,  0.1535966E+02, -0.7377693E+00 /
637
 
      data (am( 4,k,-4),k=0, 2)
638
 
     & / -0.6460216E+01, -0.4783019E+01,  0.4913297E+00 /
639
 
      data (am( 5,k,-4),k=0, 2)
640
 
     & /  0.3141830E+02, -0.3178031E+02,  0.7136013E+01 /
641
 
      data (am( 6,k,-4),k=0, 2)
642
 
     & / -0.1802509E+02,  0.1862163E+02, -0.4632843E+01 /
643
 
      data (am( 7,k,-4),k=0, 2)
644
 
     & / -0.1240412E+02,  0.2565386E+02, -0.1066570E+02 /
645
 
 
646
 
      data mexvec(-5) / 6 /
647
 
      data mlfvec(-5) / 2 /
648
 
      data ut1vec(-5) /  0.5562568E+01 /
649
 
      data ut2vec(-5) / -0.1801317E+01 /
650
 
      data alfvec(-5) /  0.4952010E-02 /
651
 
      data qmavec(-5) /  0.4500000E+01 /
652
 
      data (am( 0,k,-5),k=0, 2)
653
 
     & / -0.6031237E+01,  0.1992727E+01, -0.1076331E+01 /
654
 
      data (am( 1,k,-5),k=0, 2)
655
 
     & /  0.2933912E+01,  0.5839674E+00,  0.7509435E-01 /
656
 
      data (am( 2,k,-5),k=0, 2)
657
 
     & / -0.8284919E+01,  0.1488593E+01, -0.8251678E+00 /
658
 
      data (am( 3,k,-5),k=0, 2)
659
 
     & / -0.1925986E+02,  0.2805753E+01, -0.3015446E+01 /
660
 
      data (am( 4,k,-5),k=0, 2)
661
 
     & / -0.9480483E+01, -0.9767837E+00, -0.1165544E+01 /
662
 
      data (am( 5,k,-5),k=0, 2)
663
 
     & /  0.2193195E+02, -0.1788518E+02,  0.9460908E+01 /
664
 
      data (am( 6,k,-5),k=0, 2)
665
 
     & / -0.1327377E+02,  0.1201754E+02, -0.6277844E+01 /
666
 
 
667
 
      if(q .le. qmavec(ifl)) then
668
 
         faux5L = 0.d0
669
 
         return
670
 
      endif
671
 
 
672
 
      if(x .ge. 1.d0) then
673
 
         faux5L = 0.d0
674
 
         return
675
 
      endif
676
 
 
677
 
      tmp = log(q/alfvec(ifl))
678
 
      if(tmp .le. 0.d0) then
679
 
         faux5L = 0.d0
680
 
         return
681
 
      endif
682
 
 
683
 
      sb = log(tmp)
684
 
      sb1 = sb - 1.2d0
685
 
      sb2 = sb1*sb1
686
 
 
687
 
      do i = 0, nex
688
 
         af(i) = 0.d0
689
 
         sbx = 1.d0
690
 
         do k = 0, mlfvec(ifl)
691
 
            af(i) = af(i) + sbx*am(i,k,ifl)
692
 
            sbx = sb1*sbx
693
 
         enddo
694
 
      enddo
695
 
 
696
 
      y = -log(x)
697
 
      u = log(x/0.00001d0)
698
 
 
699
 
      part1 = af(1)*y**(1.d0+0.01d0*af(4))*(1.d0+ af(8)*u)
700
 
      part2 = af(0)*(1.d0 - x) + af(3)*x 
701
 
      part3 = x*(1.d0-x)*(af(5)+af(6)*(1.d0-x)+af(7)*x*(1.d0-x))
702
 
      part4 = ut1vec(ifl)*log(1.d0-x) + 
703
 
     &        AF(2)*log(1.d0+exp(ut2vec(ifl))-x)
704
 
 
705
 
      faux5L = exp(log(x) + part1 + part2 + part3 + part4)
706
 
 
707
 
c include threshold factor...
708
 
      faux5L = faux5L * (1.d0 - qmavec(ifl)/q)
709
 
 
710
 
      return
711
 
      end