← Back to branch summary

~maddevelopers/mg5amcnlo/3.0.2-alpha0

« back to all changes in this revision

Viewing changes to Template/Source/PDF/mrst2001.f

Added Template and HELAS into bzr

Show diffs side-by-side

added added

removed removed

Lines of Context:
Template/Source/PDF/mrst2001.f
 
1
      subroutine mrst2001(x,q,mode,upv,dnv,usea,dsea,str,chm,bot,glu)
 
2
C***************************************************************C
 
3
C                                                               C
 
4
C  This is a package for the new MRST 2001 NLO parton           C
 
5
C  distributions.                                               C     
 
6
C  Reference: A.D. Martin, R.G. Roberts, W.J. Stirling and      C
 
7
C  R.S. Thorne, hep-ph/0110215                                  C
 
8
C                                                               C
 
9
C  There are 4 pdf sets corresponding to mode = 1, 2, 3, 4      C
 
10
C                                                               C
 
11
C  Mode=1 gives the default set with Lambda(MSbar,nf=4) = 0.323 C
 
12
C  corresponding to alpha_s(M_Z) of 0.119                       C
 
13
C  This set reads a grid whose first number is 0.00927          C
 
14
C                                                               C
 
15
C  Mode=2 gives the set with Lambda(MSbar,nf=4) = 0.290         C
 
16
C  corresponding to alpha_s(M_Z) of 0.117                       C
 
17
C  This set reads a grid whose first number is 0.00953          C
 
18
C                                                               C
 
19
C  Mode=3 gives the set with Lambda(MSbar,nf=4) = 0.362         C
 
20
C  corresponding to alpha_s(M_Z) of 0.121                       C
 
21
C  This set reads a grid whose first number is 0.00889          C
 
22
C                                                               C
 
23
C  Mode=4 gives the set MRST2001J which gives better agreement  C
 
24
C  with the Tevatron inclusive jet data but has unattractive    C
 
25
C  gluon behaviour at large x (see discussion in paper)         C
 
26
C  This set has Lambda(MSbar,nf=4) = 0.353(alpha_s(M_Z) = 0.121 C 
 
27
C  This set reads a grid whose first number is 0.00826          C
 
28
C                                                               C
 
29
C   This subroutine uses an improved interpolation procedure    C 
 
30
C   for extracting values of the pdf's from the grid            C
 
31
C                                                               C
 
32
C         Comments to : W.J.Stirling@durham.ac.uk               C
 
33
C                                                               C
 
34
C***************************************************************C
 
35
      implicit real*8(a-h,o-z)
 
36
      data xmin,xmax,qsqmin,qsqmax/1d-5,1d0,1.25d0,1d7/
 
37
      q2=q*q
 
38
      if(q2.lt.qsqmin.or.q2.gt.qsqmax) print 99,q2
 
39
      if(x.lt.xmin.or.x.gt.xmax)       print 98,x
 
40
          if(mode.eq.1) then
 
41
        call mrst1(x,q2,upv,dnv,usea,dsea,str,chm,bot,glu) 
 
42
      elseif(mode.eq.2) then
 
43
        call mrst2(x,q2,upv,dnv,usea,dsea,str,chm,bot,glu) 
 
44
      elseif(mode.eq.3) then
 
45
        call mrst3(x,q2,upv,dnv,usea,dsea,str,chm,bot,glu) 
 
46
      elseif(mode.eq.4) then
 
47
        call mrst4(x,q2,upv,dnv,usea,dsea,str,chm,bot,glu)
 
48
      endif 
 
49
  99  format('  WARNING:  Q^2 VALUE IS OUT OF RANGE   ','q2= ',e10.5)
 
50
  98  format('  WARNING:   X  VALUE IS OUT OF RANGE   ','x= ',e10.5)
 
51
      return
 
52
      end
 
53
 
 
54
      subroutine mrst1(x,qsq,upv,dnv,usea,dsea,str,chm,bot,glu)
 
55
      implicit real*8(a-h,o-z)
 
56
c-fabio
 
57
      Character Tablefile*40
 
58
      data Tablefile/'alf119.dat'/
 
59
      integer IU
 
60
      common/IU/IU
 
61
c
 
62
      parameter(nx=49,nq=37,np=8,nqc0=2,nqb0=11,nqc=35,nqb=26)
 
63
      real*8 f1(nx,nq),f2(nx,nq),f3(nx,nq),f4(nx,nq),f5(nx,nq),
 
64
     .f6(nx,nq),f7(nx,nq),f8(nx,nq),fc(nx,nqc),fb(nx,nqb)
 
65
      real*8 qq(nq),xx(nx),cc1(nx,nq,4,4),cc2(nx,nq,4,4),
 
66
     .cc3(nx,nq,4,4),cc4(nx,nq,4,4),cc6(nx,nq,4,4),cc8(nx,nq,4,4),
 
67
     .ccc(nx,nqc,4,4),ccb(nx,nqb,4,4)
 
68
      real*8 xxl(nx),qql(nq),qqlc(nqc),qqlb(nqb)
 
69
      data xx/1d-5,2d-5,4d-5,6d-5,8d-5,
 
70
     .        1d-4,2d-4,4d-4,6d-4,8d-4,
 
71
     .        1d-3,2d-3,4d-3,6d-3,8d-3,
 
72
     .        1d-2,1.4d-2,2d-2,3d-2,4d-2,6d-2,8d-2,
 
73
     .     .1d0,.125d0,.15d0,.175d0,.2d0,.225d0,.25d0,.275d0,
 
74
     .     .3d0,.325d0,.35d0,.375d0,.4d0,.425d0,.45d0,.475d0,
 
75
     .     .5d0,.525d0,.55d0,.575d0,.6d0,.65d0,.7d0,.75d0,
 
76
     .     .8d0,.9d0,1d0/
 
77
      data qq/1.25d0,1.5d0,2d0,2.5d0,3.2d0,4d0,5d0,6.4d0,8d0,1d1,
 
78
     .        1.2d1,1.8d1,2.6d1,4d1,6.4d1,1d2,
 
79
     .        1.6d2,2.4d2,4d2,6.4d2,1d3,1.8d3,3.2d3,5.6d3,1d4,
 
80
     .        1.8d4,3.2d4,5.6d4,1d5,1.8d5,3.2d5,5.6d5,1d6,
 
81
     .        1.8d6,3.2d6,5.6d6,1d7/
 
82
c      data xmin,xmax,qsqmin,qsqmax/1d-5,1d0,1.25d0,1d7/
 
83
      data init/0/
 
84
      save
 
85
      xsave=x
 
86
      q2save=qsq
 
87
      if(init.ne.0) goto 10
 
88
c        open(unit=33,file='Pdfdata/alf119.dat',status='old')
 
89
        call OpenData(Tablefile)
 
90
        do 20 n=1,nx-1
 
91
        do 20 m=1,nq
 
92
        read(IU,50)f1(n,m),f2(n,m),f3(n,m),f4(n,m),
 
93
     .            f5(n,m),f7(n,m),f6(n,m),f8(n,m)
 
94
c notation: 1=uval 2=val 3=glue 4=usea 5=chm 6=str 7=btm 8=dsea
 
95
  20  continue
 
96
      do 40 m=1,nq
 
97
      f1(nx,m)=0.d0
 
98
      f2(nx,m)=0.d0
 
99
      f3(nx,m)=0.d0
 
100
      f4(nx,m)=0.d0
 
101
      f5(nx,m)=0.d0
 
102
      f6(nx,m)=0.d0
 
103
      f7(nx,m)=0.d0
 
104
      f8(nx,m)=0.d0
 
105
  40  continue
 
106
      do n=1,nx
 
107
      xxl(n)=dlog(xx(n))
 
108
      enddo
 
109
      do m=1,nq
 
110
      qql(m)=dlog(qq(m))
 
111
      enddo
 
112
 
 
113
      call jeppe1(nx,nq,xxl,qql,f1,cc1)
 
114
      call jeppe1(nx,nq,xxl,qql,f2,cc2)
 
115
      call jeppe1(nx,nq,xxl,qql,f3,cc3)
 
116
      call jeppe1(nx,nq,xxl,qql,f4,cc4)
 
117
      call jeppe1(nx,nq,xxl,qql,f6,cc6)
 
118
      call jeppe1(nx,nq,xxl,qql,f8,cc8)
 
119
 
 
120
      emc2=2.045d0
 
121
      emb2=18.5d0
 
122
 
 
123
      do 44 m=1,nqc
 
124
      qqlc(m)=qql(m+nqc0)
 
125
      do 44 n=1,nx
 
126
      fc(n,m)=f5(n,m+nqc0)
 
127
   44 continue
 
128
      qqlc(1)=dlog(emc2)
 
129
      call jeppe1(nx,nqc,xxl,qqlc,fc,ccc)
 
130
 
 
131
      do 45 m=1,nqb
 
132
      qqlb(m)=qql(m+nqb0)
 
133
      do 45 n=1,nx
 
134
      fb(n,m)=f7(n,m+nqb0)
 
135
   45 continue
 
136
      qqlb(1)=dlog(emb2)
 
137
      call jeppe1(nx,nqb,xxl,qqlb,fb,ccb)
 
138
 
 
139
 
 
140
      init=1
 
141
   10 continue
 
142
      
 
143
      xlog=dlog(x)
 
144
      qsqlog=dlog(qsq)
 
145
 
 
146
      call jeppe2(xlog,qsqlog,nx,nq,xxl,qql,cc1,upv)
 
147
      call jeppe2(xlog,qsqlog,nx,nq,xxl,qql,cc2,dnv)
 
148
      call jeppe2(xlog,qsqlog,nx,nq,xxl,qql,cc3,glu)
 
149
      call jeppe2(xlog,qsqlog,nx,nq,xxl,qql,cc4,usea)
 
150
      call jeppe2(xlog,qsqlog,nx,nq,xxl,qql,cc6,str)
 
151
      call jeppe2(xlog,qsqlog,nx,nq,xxl,qql,cc8,dsea)
 
152
 
 
153
      chm=0.d0
 
154
      if(qsq.gt.emc2) then 
 
155
      call jeppe2(xlog,qsqlog,nx,nqc,xxl,qqlc,ccc,chm)
 
156
      endif
 
157
 
 
158
      bot=0.d0
 
159
      if(qsq.gt.emb2) then 
 
160
      call jeppe2(xlog,qsqlog,nx,nqb,xxl,qqlb,ccb,bot)
 
161
      endif
 
162
 
 
163
      x=xsave
 
164
      qsq=q2save
 
165
      return
 
166
   50 format(8f10.5)
 
167
      end
 
168
 
 
169
      subroutine mrst2(x,qsq,upv,dnv,usea,dsea,str,chm,bot,glu)
 
170
      implicit real*8(a-h,o-z)
 
171
c-fabio
 
172
      Character Tablefile*40
 
173
      data Tablefile/'alf117.dat'/
 
174
      integer IU
 
175
      common/IU/IU
 
176
c
 
177
      parameter(nx=49,nq=37,np=8,nqc0=2,nqb0=11,nqc=35,nqb=26)
 
178
      real*8 f1(nx,nq),f2(nx,nq),f3(nx,nq),f4(nx,nq),f5(nx,nq),
 
179
     .f6(nx,nq),f7(nx,nq),f8(nx,nq),fc(nx,nqc),fb(nx,nqb)
 
180
      real*8 qq(nq),xx(nx),cc1(nx,nq,4,4),cc2(nx,nq,4,4),
 
181
     .cc3(nx,nq,4,4),cc4(nx,nq,4,4),cc6(nx,nq,4,4),cc8(nx,nq,4,4),
 
182
     .ccc(nx,nqc,4,4),ccb(nx,nqb,4,4)
 
183
      real*8 xxl(nx),qql(nq),qqlc(nqc),qqlb(nqb)
 
184
      data xx/1d-5,2d-5,4d-5,6d-5,8d-5,
 
185
     .        1d-4,2d-4,4d-4,6d-4,8d-4,
 
186
     .        1d-3,2d-3,4d-3,6d-3,8d-3,
 
187
     .        1d-2,1.4d-2,2d-2,3d-2,4d-2,6d-2,8d-2,
 
188
     .     .1d0,.125d0,.15d0,.175d0,.2d0,.225d0,.25d0,.275d0,
 
189
     .     .3d0,.325d0,.35d0,.375d0,.4d0,.425d0,.45d0,.475d0,
 
190
     .     .5d0,.525d0,.55d0,.575d0,.6d0,.65d0,.7d0,.75d0,
 
191
     .     .8d0,.9d0,1d0/
 
192
      data qq/1.25d0,1.5d0,2d0,2.5d0,3.2d0,4d0,5d0,6.4d0,8d0,1d1,
 
193
     .        1.2d1,1.8d1,2.6d1,4d1,6.4d1,1d2,
 
194
     .        1.6d2,2.4d2,4d2,6.4d2,1d3,1.8d3,3.2d3,5.6d3,1d4,
 
195
     .        1.8d4,3.2d4,5.6d4,1d5,1.8d5,3.2d5,5.6d5,1d6,
 
196
     .        1.8d6,3.2d6,5.6d6,1d7/
 
197
c      data xmin,xmax,qsqmin,qsqmax/1d-5,1d0,1.25d0,1d7/
 
198
      data init/0/
 
199
      save
 
200
      xsave=x
 
201
      q2save=qsq
 
202
      if(init.ne.0) goto 10
 
203
c        open(unit=33,file='Pdfdata/alf117.dat',status='old')
 
204
        call OpenData(Tablefile)
 
205
        do 20 n=1,nx-1
 
206
        do 20 m=1,nq
 
207
        read(IU,50)f1(n,m),f2(n,m),f3(n,m),f4(n,m),
 
208
     .            f5(n,m),f7(n,m),f6(n,m),f8(n,m)
 
209
c notation: 1=uval 2=val 3=glue 4=usea 5=chm 6=str 7=btm 8=dsea
 
210
  20  continue
 
211
      do 40 m=1,nq
 
212
      f1(nx,m)=0.d0
 
213
      f2(nx,m)=0.d0
 
214
      f3(nx,m)=0.d0
 
215
      f4(nx,m)=0.d0
 
216
      f5(nx,m)=0.d0
 
217
      f6(nx,m)=0.d0
 
218
      f7(nx,m)=0.d0
 
219
      f8(nx,m)=0.d0
 
220
  40  continue
 
221
      do n=1,nx
 
222
      xxl(n)=dlog(xx(n))
 
223
      enddo
 
224
      do m=1,nq
 
225
      qql(m)=dlog(qq(m))
 
226
      enddo
 
227
 
 
228
      call jeppe1(nx,nq,xxl,qql,f1,cc1)
 
229
      call jeppe1(nx,nq,xxl,qql,f2,cc2)
 
230
      call jeppe1(nx,nq,xxl,qql,f3,cc3)
 
231
      call jeppe1(nx,nq,xxl,qql,f4,cc4)
 
232
      call jeppe1(nx,nq,xxl,qql,f6,cc6)
 
233
      call jeppe1(nx,nq,xxl,qql,f8,cc8)
 
234
 
 
235
      emc2=2.045d0
 
236
      emb2=18.5d0
 
237
 
 
238
      do 44 m=1,nqc
 
239
      qqlc(m)=qql(m+nqc0)
 
240
      do 44 n=1,nx
 
241
      fc(n,m)=f5(n,m+nqc0)
 
242
   44 continue
 
243
      qqlc(1)=dlog(emc2)
 
244
      call jeppe1(nx,nqc,xxl,qqlc,fc,ccc)
 
245
 
 
246
      do 45 m=1,nqb
 
247
      qqlb(m)=qql(m+nqb0)
 
248
      do 45 n=1,nx
 
249
      fb(n,m)=f7(n,m+nqb0)
 
250
   45 continue
 
251
      qqlb(1)=dlog(emb2)
 
252
      call jeppe1(nx,nqb,xxl,qqlb,fb,ccb)
 
253
 
 
254
 
 
255
      init=1
 
256
   10 continue
 
257
      
 
258
      xlog=dlog(x)
 
259
      qsqlog=dlog(qsq)
 
260
 
 
261
      call jeppe2(xlog,qsqlog,nx,nq,xxl,qql,cc1,upv)
 
262
      call jeppe2(xlog,qsqlog,nx,nq,xxl,qql,cc2,dnv)
 
263
      call jeppe2(xlog,qsqlog,nx,nq,xxl,qql,cc3,glu)
 
264
      call jeppe2(xlog,qsqlog,nx,nq,xxl,qql,cc4,usea)
 
265
      call jeppe2(xlog,qsqlog,nx,nq,xxl,qql,cc6,str)
 
266
      call jeppe2(xlog,qsqlog,nx,nq,xxl,qql,cc8,dsea)
 
267
 
 
268
      chm=0.d0
 
269
      if(qsq.gt.emc2) then 
 
270
      call jeppe2(xlog,qsqlog,nx,nqc,xxl,qqlc,ccc,chm)
 
271
      endif
 
272
 
 
273
      bot=0.d0
 
274
      if(qsq.gt.emb2) then 
 
275
      call jeppe2(xlog,qsqlog,nx,nqb,xxl,qqlb,ccb,bot)
 
276
      endif
 
277
 
 
278
      x=xsave
 
279
      qsq=q2save
 
280
      return
 
281
   50 format(8f10.5)
 
282
      end
 
283
 
 
284
      subroutine mrst3(x,qsq,upv,dnv,usea,dsea,str,chm,bot,glu)
 
285
      implicit real*8(a-h,o-z)
 
286
c-fabio
 
287
      Character Tablefile*40
 
288
      data Tablefile/'alf121.dat'/
 
289
      integer IU
 
290
      common/IU/IU
 
291
c
 
292
      parameter(nx=49,nq=37,np=8,nqc0=2,nqb0=11,nqc=35,nqb=26)
 
293
      real*8 f1(nx,nq),f2(nx,nq),f3(nx,nq),f4(nx,nq),f5(nx,nq),
 
294
     .f6(nx,nq),f7(nx,nq),f8(nx,nq),fc(nx,nqc),fb(nx,nqb)
 
295
      real*8 qq(nq),xx(nx),cc1(nx,nq,4,4),cc2(nx,nq,4,4),
 
296
     .cc3(nx,nq,4,4),cc4(nx,nq,4,4),cc6(nx,nq,4,4),cc8(nx,nq,4,4),
 
297
     .ccc(nx,nqc,4,4),ccb(nx,nqb,4,4)
 
298
      real*8 xxl(nx),qql(nq),qqlc(nqc),qqlb(nqb)
 
299
      data xx/1d-5,2d-5,4d-5,6d-5,8d-5,
 
300
     .        1d-4,2d-4,4d-4,6d-4,8d-4,
 
301
     .        1d-3,2d-3,4d-3,6d-3,8d-3,
 
302
     .        1d-2,1.4d-2,2d-2,3d-2,4d-2,6d-2,8d-2,
 
303
     .     .1d0,.125d0,.15d0,.175d0,.2d0,.225d0,.25d0,.275d0,
 
304
     .     .3d0,.325d0,.35d0,.375d0,.4d0,.425d0,.45d0,.475d0,
 
305
     .     .5d0,.525d0,.55d0,.575d0,.6d0,.65d0,.7d0,.75d0,
 
306
     .     .8d0,.9d0,1d0/
 
307
      data qq/1.25d0,1.5d0,2d0,2.5d0,3.2d0,4d0,5d0,6.4d0,8d0,1d1,
 
308
     .        1.2d1,1.8d1,2.6d1,4d1,6.4d1,1d2,
 
309
     .        1.6d2,2.4d2,4d2,6.4d2,1d3,1.8d3,3.2d3,5.6d3,1d4,
 
310
     .        1.8d4,3.2d4,5.6d4,1d5,1.8d5,3.2d5,5.6d5,1d6,
 
311
     .        1.8d6,3.2d6,5.6d6,1d7/
 
312
c      data xmin,xmax,qsqmin,qsqmax/1d-5,1d0,1.25d0,1d7/
 
313
      data init/0/
 
314
      save
 
315
      xsave=x
 
316
      q2save=qsq
 
317
      if(init.ne.0) goto 10
 
318
c        open(unit=33,file='Pdfdata/alf121.dat',status='old')
 
319
        call OpenData(Tablefile)
 
320
        do 20 n=1,nx-1
 
321
        do 20 m=1,nq
 
322
        read(IU,50)f1(n,m),f2(n,m),f3(n,m),f4(n,m),
 
323
     .            f5(n,m),f7(n,m),f6(n,m),f8(n,m)
 
324
c notation: 1=uval 2=val 3=glue 4=usea 5=chm 6=str 7=btm 8=dsea
 
325
  20  continue
 
326
      do 40 m=1,nq
 
327
      f1(nx,m)=0.d0
 
328
      f2(nx,m)=0.d0
 
329
      f3(nx,m)=0.d0
 
330
      f4(nx,m)=0.d0
 
331
      f5(nx,m)=0.d0
 
332
      f6(nx,m)=0.d0
 
333
      f7(nx,m)=0.d0
 
334
      f8(nx,m)=0.d0
 
335
  40  continue
 
336
      do n=1,nx
 
337
      xxl(n)=dlog(xx(n))
 
338
      enddo
 
339
      do m=1,nq
 
340
      qql(m)=dlog(qq(m))
 
341
      enddo
 
342
 
 
343
      call jeppe1(nx,nq,xxl,qql,f1,cc1)
 
344
      call jeppe1(nx,nq,xxl,qql,f2,cc2)
 
345
      call jeppe1(nx,nq,xxl,qql,f3,cc3)
 
346
      call jeppe1(nx,nq,xxl,qql,f4,cc4)
 
347
      call jeppe1(nx,nq,xxl,qql,f6,cc6)
 
348
      call jeppe1(nx,nq,xxl,qql,f8,cc8)
 
349
 
 
350
      emc2=2.045d0
 
351
      emb2=18.5d0
 
352
 
 
353
      do 44 m=1,nqc
 
354
      qqlc(m)=qql(m+nqc0)
 
355
      do 44 n=1,nx
 
356
      fc(n,m)=f5(n,m+nqc0)
 
357
   44 continue
 
358
      qqlc(1)=dlog(emc2)
 
359
      call jeppe1(nx,nqc,xxl,qqlc,fc,ccc)
 
360
 
 
361
      do 45 m=1,nqb
 
362
      qqlb(m)=qql(m+nqb0)
 
363
      do 45 n=1,nx
 
364
      fb(n,m)=f7(n,m+nqb0)
 
365
   45 continue
 
366
      qqlb(1)=dlog(emb2)
 
367
      call jeppe1(nx,nqb,xxl,qqlb,fb,ccb)
 
368
 
 
369
 
 
370
      init=1
 
371
   10 continue
 
372
      
 
373
      xlog=dlog(x)
 
374
      qsqlog=dlog(qsq)
 
375
 
 
376
      call jeppe2(xlog,qsqlog,nx,nq,xxl,qql,cc1,upv)
 
377
      call jeppe2(xlog,qsqlog,nx,nq,xxl,qql,cc2,dnv)
 
378
      call jeppe2(xlog,qsqlog,nx,nq,xxl,qql,cc3,glu)
 
379
      call jeppe2(xlog,qsqlog,nx,nq,xxl,qql,cc4,usea)
 
380
      call jeppe2(xlog,qsqlog,nx,nq,xxl,qql,cc6,str)
 
381
      call jeppe2(xlog,qsqlog,nx,nq,xxl,qql,cc8,dsea)
 
382
 
 
383
      chm=0.d0
 
384
      if(qsq.gt.emc2) then 
 
385
      call jeppe2(xlog,qsqlog,nx,nqc,xxl,qqlc,ccc,chm)
 
386
      endif
 
387
 
 
388
      bot=0.d0
 
389
      if(qsq.gt.emb2) then 
 
390
      call jeppe2(xlog,qsqlog,nx,nqb,xxl,qqlb,ccb,bot)
 
391
      endif
 
392
 
 
393
      x=xsave
 
394
      qsq=q2save
 
395
      return
 
396
   50 format(8f10.5)
 
397
      end
 
398
 
 
399
      subroutine mrst4(x,qsq,upv,dnv,usea,dsea,str,chm,bot,glu)
 
400
      implicit real*8(a-h,o-z)
 
401
c-fabio
 
402
      Character Tablefile*40
 
403
      data Tablefile/'j121.dat'/
 
404
      integer IU
 
405
      common/IU/IU
 
406
c
 
407
      parameter(nx=49,nq=37,np=8,nqc0=2,nqb0=11,nqc=35,nqb=26)
 
408
      real*8 f1(nx,nq),f2(nx,nq),f3(nx,nq),f4(nx,nq),f5(nx,nq),
 
409
     .f6(nx,nq),f7(nx,nq),f8(nx,nq),fc(nx,nqc),fb(nx,nqb)
 
410
      real*8 qq(nq),xx(nx),cc1(nx,nq,4,4),cc2(nx,nq,4,4),
 
411
     .cc3(nx,nq,4,4),cc4(nx,nq,4,4),cc6(nx,nq,4,4),cc8(nx,nq,4,4),
 
412
     .ccc(nx,nqc,4,4),ccb(nx,nqb,4,4)
 
413
      real*8 xxl(nx),qql(nq),qqlc(nqc),qqlb(nqb)
 
414
      data xx/1d-5,2d-5,4d-5,6d-5,8d-5,
 
415
     .        1d-4,2d-4,4d-4,6d-4,8d-4,
 
416
     .        1d-3,2d-3,4d-3,6d-3,8d-3,
 
417
     .        1d-2,1.4d-2,2d-2,3d-2,4d-2,6d-2,8d-2,
 
418
     .     .1d0,.125d0,.15d0,.175d0,.2d0,.225d0,.25d0,.275d0,
 
419
     .     .3d0,.325d0,.35d0,.375d0,.4d0,.425d0,.45d0,.475d0,
 
420
     .     .5d0,.525d0,.55d0,.575d0,.6d0,.65d0,.7d0,.75d0,
 
421
     .     .8d0,.9d0,1d0/
 
422
      data qq/1.25d0,1.5d0,2d0,2.5d0,3.2d0,4d0,5d0,6.4d0,8d0,1d1,
 
423
     .        1.2d1,1.8d1,2.6d1,4d1,6.4d1,1d2,
 
424
     .        1.6d2,2.4d2,4d2,6.4d2,1d3,1.8d3,3.2d3,5.6d3,1d4,
 
425
     .        1.8d4,3.2d4,5.6d4,1d5,1.8d5,3.2d5,5.6d5,1d6,
 
426
     .        1.8d6,3.2d6,5.6d6,1d7/
 
427
c      data xmin,xmax,qsqmin,qsqmax/1d-5,1d0,1.25d0,1d7/
 
428
      data init/0/
 
429
      save
 
430
      xsave=x
 
431
      q2save=qsq
 
432
      if(init.ne.0) goto 10
 
433
c        open(unit=33,file='Pdfdata/j121.dat',status='old')
 
434
        call OpenData(Tablefile)
 
435
        do 20 n=1,nx-1
 
436
        do 20 m=1,nq
 
437
        read(IU,50)f1(n,m),f2(n,m),f3(n,m),f4(n,m),
 
438
     .            f5(n,m),f7(n,m),f6(n,m),f8(n,m)
 
439
c notation: 1=uval 2=val 3=glue 4=usea 5=chm 6=str 7=btm 8=dsea
 
440
  20  continue
 
441
      do 40 m=1,nq
 
442
      f1(nx,m)=0.d0
 
443
      f2(nx,m)=0.d0
 
444
      f3(nx,m)=0.d0
 
445
      f4(nx,m)=0.d0
 
446
      f5(nx,m)=0.d0
 
447
      f6(nx,m)=0.d0
 
448
      f7(nx,m)=0.d0
 
449
      f8(nx,m)=0.d0
 
450
  40  continue
 
451
      do n=1,nx
 
452
      xxl(n)=dlog(xx(n))
 
453
      enddo
 
454
      do m=1,nq
 
455
      qql(m)=dlog(qq(m))
 
456
      enddo
 
457
 
 
458
      call jeppe1(nx,nq,xxl,qql,f1,cc1)
 
459
      call jeppe1(nx,nq,xxl,qql,f2,cc2)
 
460
      call jeppe1(nx,nq,xxl,qql,f3,cc3)
 
461
      call jeppe1(nx,nq,xxl,qql,f4,cc4)
 
462
      call jeppe1(nx,nq,xxl,qql,f6,cc6)
 
463
      call jeppe1(nx,nq,xxl,qql,f8,cc8)
 
464
 
 
465
      emc2=2.045d0
 
466
      emb2=18.5d0
 
467
 
 
468
      do 44 m=1,nqc
 
469
      qqlc(m)=qql(m+nqc0)
 
470
      do 44 n=1,nx
 
471
      fc(n,m)=f5(n,m+nqc0)
 
472
   44 continue
 
473
      qqlc(1)=dlog(emc2)
 
474
      call jeppe1(nx,nqc,xxl,qqlc,fc,ccc)
 
475
 
 
476
      do 45 m=1,nqb
 
477
      qqlb(m)=qql(m+nqb0)
 
478
      do 45 n=1,nx
 
479
      fb(n,m)=f7(n,m+nqb0)
 
480
   45 continue
 
481
      qqlb(1)=dlog(emb2)
 
482
      call jeppe1(nx,nqb,xxl,qqlb,fb,ccb)
 
483
 
 
484
 
 
485
      init=1
 
486
   10 continue
 
487
      
 
488
      xlog=dlog(x)
 
489
      qsqlog=dlog(qsq)
 
490
 
 
491
      call jeppe2(xlog,qsqlog,nx,nq,xxl,qql,cc1,upv)
 
492
      call jeppe2(xlog,qsqlog,nx,nq,xxl,qql,cc2,dnv)
 
493
      call jeppe2(xlog,qsqlog,nx,nq,xxl,qql,cc3,glu)
 
494
      call jeppe2(xlog,qsqlog,nx,nq,xxl,qql,cc4,usea)
 
495
      call jeppe2(xlog,qsqlog,nx,nq,xxl,qql,cc6,str)
 
496
      call jeppe2(xlog,qsqlog,nx,nq,xxl,qql,cc8,dsea)
 
497
 
 
498
      chm=0.d0
 
499
      if(qsq.gt.emc2) then 
 
500
      call jeppe2(xlog,qsqlog,nx,nqc,xxl,qqlc,ccc,chm)
 
501
      endif
 
502
 
 
503
      bot=0.d0
 
504
      if(qsq.gt.emb2) then 
 
505
      call jeppe2(xlog,qsqlog,nx,nqb,xxl,qqlb,ccb,bot)
 
506
      endif
 
507
 
 
508
      x=xsave
 
509
      qsq=q2save
 
510
      return
 
511
   50 format(8f10.5)
 
512
      end
 
513