← Back to branch summary

~madteam/mg5amcnlo/series2.0

« back to all changes in this revision

Viewing changes to Template/FKS-born/Source/PDF/mrst2002.f

  • Committer: Marco Zaro
  • Date: 2012-08-28 21:06:34 UTC
  • mto: (78.35.14 AutoMint)
  • mto: This revision was merged to the branch mainline in revision 249.
  • Revision ID: marco.zaro@gmail.com-20120828210634-5a06yvda3hplw8ur
doing some renaming:
 Template/FKS-born -> Template/NLO
 fks_born.py -> fks_base.py
 fks_born_helas_objects.py -> fks_helas_objects.py
 export_fks_born.py -> export_fks.py

also functions/classes and tests renamed
all unit tests ok, exporting ok

Show diffs side-by-side

added added

removed removed

Lines of Context:
Template/FKS-born/Source/PDF/mrst2002.f
1
 
      subroutine mrst2002(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 2002 updated NLO and      C
5
 
C  NNLO parton distributions.                                   C 
6
 
C  Reference: A.D. Martin, R.G. Roberts, W.J. Stirling and      C
7
 
C  R.S. Thorne, hep-ph/0211080                                  C
8
 
C                                                               C
9
 
C  There are 2 pdf sets corresponding to mode = 1, 2            C
10
 
C                                                               C
11
 
C  Mode=1 gives the NLO set with alpha_s(M_Z,NLO) = 0.1197      C  
12
 
C  This set reads a grid whose first number is 0.00949          C
13
 
C                                                               C
14
 
C  Mode=2 gives the NNLO set with alpha_s(M_Z,NNLO) = 0.1154    C
15
 
C  This set reads a grid whose first number is 0.00685          C
16
 
C                                                               C
17
 
C         Comments to : W.J.Stirling@durham.ac.uk               C
18
 
C                                                               C
19
 
C***************************************************************C
20
 
      implicit real*8(a-h,o-z)
21
 
      data xmin,xmax,qsqmin,qsqmax/1d-5,1d0,1.25d0,1d7/
22
 
      q2=q*q
23
 
      if(q2.lt.qsqmin.or.q2.gt.qsqmax) print 99,q2
24
 
      if(x.lt.xmin.or.x.gt.xmax)       print 98,x
25
 
          if(mode.eq.1) then
26
 
        call mrst_02_1(x,q2,upv,dnv,usea,dsea,str,chm,bot,glu) 
27
 
      elseif(mode.eq.2) then
28
 
        call mrst_02_1(x,q2,upv,dnv,usea,dsea,str,chm,bot,glu) 
29
 
      endif 
30
 
  99  format('  WARNING:  Q^2 VALUE IS OUT OF RANGE   ','q2= ',e10.5)
31
 
  98  format('  WARNING:   X  VALUE IS OUT OF RANGE   ','x= ',e10.5)
32
 
      return
33
 
      end
34
 
 
35
 
      subroutine mrst_02_1(x,qsq,upv,dnv,usea,dsea,str,chm,bot,glu)
36
 
      implicit real*8(a-h,o-z)
37
 
c-fabio
38
 
      Character Tablefile*40
39
 
      data Tablefile/'mrst2002nlo.dat'/
40
 
      integer IU
41
 
      common/IU/IU
42
 
c
43
 
      parameter(nx=49,nq=37,np=8,nqc0=2,nqb0=11,nqc=35,nqb=26)
44
 
      real*8 f1(nx,nq),f2(nx,nq),f3(nx,nq),f4(nx,nq),f5(nx,nq),
45
 
     .f6(nx,nq),f7(nx,nq),f8(nx,nq),fc(nx,nqc),fb(nx,nqb)
46
 
      real*8 qq(nq),xx(nx),cc1(nx,nq,4,4),cc2(nx,nq,4,4),
47
 
     .cc3(nx,nq,4,4),cc4(nx,nq,4,4),cc6(nx,nq,4,4),cc8(nx,nq,4,4),
48
 
     .ccc(nx,nqc,4,4),ccb(nx,nqb,4,4)
49
 
      real*8 xxl(nx),qql(nq),qqlc(nqc),qqlb(nqb)
50
 
      data xx/1d-5,2d-5,4d-5,6d-5,8d-5,
51
 
     .        1d-4,2d-4,4d-4,6d-4,8d-4,
52
 
     .        1d-3,2d-3,4d-3,6d-3,8d-3,
53
 
     .        1d-2,1.4d-2,2d-2,3d-2,4d-2,6d-2,8d-2,
54
 
     .     .1d0,.125d0,.15d0,.175d0,.2d0,.225d0,.25d0,.275d0,
55
 
     .     .3d0,.325d0,.35d0,.375d0,.4d0,.425d0,.45d0,.475d0,
56
 
     .     .5d0,.525d0,.55d0,.575d0,.6d0,.65d0,.7d0,.75d0,
57
 
     .     .8d0,.9d0,1d0/
58
 
      data qq/1.25d0,1.5d0,2d0,2.5d0,3.2d0,4d0,5d0,6.4d0,8d0,1d1,
59
 
     .        1.2d1,1.8d1,2.6d1,4d1,6.4d1,1d2,
60
 
     .        1.6d2,2.4d2,4d2,6.4d2,1d3,1.8d3,3.2d3,5.6d3,1d4,
61
 
     .        1.8d4,3.2d4,5.6d4,1d5,1.8d5,3.2d5,5.6d5,1d6,
62
 
     .        1.8d6,3.2d6,5.6d6,1d7/
63
 
      data xmin,xmax,qsqmin,qsqmax/1d-5,1d0,1.25d0,1d7/
64
 
      data init/0/
65
 
      save
66
 
      xsave=x
67
 
      q2save=qsq
68
 
      if(init.ne.0) goto 10
69
 
c        open(unit=33,file='mrst2002nlo.dat',status='old')
70
 
        call OpenData(Tablefile)
71
 
        do 20 n=1,nx-1
72
 
        do 20 m=1,nq
73
 
        read(IU,50)f1(n,m),f2(n,m),f3(n,m),f4(n,m),
74
 
     .            f5(n,m),f7(n,m),f6(n,m),f8(n,m)
75
 
c notation: 1=uval 2=val 3=glue 4=usea 5=chm 6=str 7=btm 8=dsea
76
 
  20  continue
77
 
      do 40 m=1,nq
78
 
      f1(nx,m)=0.d0
79
 
      f2(nx,m)=0.d0
80
 
      f3(nx,m)=0.d0
81
 
      f4(nx,m)=0.d0
82
 
      f5(nx,m)=0.d0
83
 
      f6(nx,m)=0.d0
84
 
      f7(nx,m)=0.d0
85
 
      f8(nx,m)=0.d0
86
 
  40  continue
87
 
      do n=1,nx
88
 
      xxl(n)=dlog(xx(n))
89
 
      enddo
90
 
      do m=1,nq
91
 
      qql(m)=dlog(qq(m))
92
 
      enddo
93
 
 
94
 
      call jeppe1(nx,nq,xxl,qql,f1,cc1)
95
 
      call jeppe1(nx,nq,xxl,qql,f2,cc2)
96
 
      call jeppe1(nx,nq,xxl,qql,f3,cc3)
97
 
      call jeppe1(nx,nq,xxl,qql,f4,cc4)
98
 
      call jeppe1(nx,nq,xxl,qql,f6,cc6)
99
 
      call jeppe1(nx,nq,xxl,qql,f8,cc8)
100
 
 
101
 
      emc2=2.045
102
 
      emb2=18.5
103
 
 
104
 
      do 44 m=1,nqc
105
 
      qqlc(m)=qql(m+nqc0)
106
 
      do 44 n=1,nx
107
 
      fc(n,m)=f5(n,m+nqc0)
108
 
   44 continue
109
 
      qqlc(1)=dlog(emc2)
110
 
      call jeppe1(nx,nqc,xxl,qqlc,fc,ccc)
111
 
 
112
 
      do 45 m=1,nqb
113
 
      qqlb(m)=qql(m+nqb0)
114
 
      do 45 n=1,nx
115
 
      fb(n,m)=f7(n,m+nqb0)
116
 
   45 continue
117
 
      qqlb(1)=dlog(emb2)
118
 
      call jeppe1(nx,nqb,xxl,qqlb,fb,ccb)
119
 
 
120
 
 
121
 
      init=1
122
 
   10 continue
123
 
      
124
 
      xlog=dlog(x)
125
 
      qsqlog=dlog(qsq)
126
 
 
127
 
      call jeppe2(xlog,qsqlog,nx,nq,xxl,qql,cc1,upv)
128
 
      call jeppe2(xlog,qsqlog,nx,nq,xxl,qql,cc2,dnv)
129
 
      call jeppe2(xlog,qsqlog,nx,nq,xxl,qql,cc3,glu)
130
 
      call jeppe2(xlog,qsqlog,nx,nq,xxl,qql,cc4,usea)
131
 
      call jeppe2(xlog,qsqlog,nx,nq,xxl,qql,cc6,str)
132
 
      call jeppe2(xlog,qsqlog,nx,nq,xxl,qql,cc8,dsea)
133
 
 
134
 
      chm=0.d0
135
 
      if(qsq.gt.emc2) then 
136
 
      call jeppe2(xlog,qsqlog,nx,nqc,xxl,qqlc,ccc,chm)
137
 
      endif
138
 
 
139
 
      bot=0.d0
140
 
      if(qsq.gt.emb2) then 
141
 
      call jeppe2(xlog,qsqlog,nx,nqb,xxl,qqlb,ccb,bot)
142
 
      endif
143
 
 
144
 
      x=xsave
145
 
      qsq=q2save
146
 
      return
147
 
   50 format(8f10.5)
148
 
      end
149
 
 
150
 
      subroutine mrst_02_2(x,qsq,upv,dnv,usea,dsea,str,chm,bot,glu)
151
 
      implicit real*8(a-h,o-z)
152
 
c-fabio
153
 
      Character Tablefile*40
154
 
      data Tablefile/'mrst2002nnlo.dat'/
155
 
      integer IU
156
 
      common/IU/IU
157
 
c
158
 
      parameter(nx=49,nq=37,np=8,nqc0=2,nqb0=11,nqc=35,nqb=26)
159
 
      real*8 f1(nx,nq),f2(nx,nq),f3(nx,nq),f4(nx,nq),f5(nx,nq),
160
 
     .f6(nx,nq),f7(nx,nq),f8(nx,nq),fc(nx,nqc),fb(nx,nqb)
161
 
      real*8 qq(nq),xx(nx),cc1(nx,nq,4,4),cc2(nx,nq,4,4),
162
 
     .cc3(nx,nq,4,4),cc4(nx,nq,4,4),cc6(nx,nq,4,4),cc8(nx,nq,4,4),
163
 
     .ccc(nx,nqc,4,4),ccb(nx,nqb,4,4)
164
 
      real*8 xxl(nx),qql(nq),qqlc(nqc),qqlb(nqb)
165
 
      data xx/1d-5,2d-5,4d-5,6d-5,8d-5,
166
 
     .        1d-4,2d-4,4d-4,6d-4,8d-4,
167
 
     .        1d-3,2d-3,4d-3,6d-3,8d-3,
168
 
     .        1d-2,1.4d-2,2d-2,3d-2,4d-2,6d-2,8d-2,
169
 
     .     .1d0,.125d0,.15d0,.175d0,.2d0,.225d0,.25d0,.275d0,
170
 
     .     .3d0,.325d0,.35d0,.375d0,.4d0,.425d0,.45d0,.475d0,
171
 
     .     .5d0,.525d0,.55d0,.575d0,.6d0,.65d0,.7d0,.75d0,
172
 
     .     .8d0,.9d0,1d0/
173
 
      data qq/1.25d0,1.5d0,2d0,2.5d0,3.2d0,4d0,5d0,6.4d0,8d0,1d1,
174
 
     .        1.2d1,1.8d1,2.6d1,4d1,6.4d1,1d2,
175
 
     .        1.6d2,2.4d2,4d2,6.4d2,1d3,1.8d3,3.2d3,5.6d3,1d4,
176
 
     .        1.8d4,3.2d4,5.6d4,1d5,1.8d5,3.2d5,5.6d5,1d6,
177
 
     .        1.8d6,3.2d6,5.6d6,1d7/
178
 
      data xmin,xmax,qsqmin,qsqmax/1d-5,1d0,1.25d0,1d7/
179
 
      data init/0/
180
 
      save
181
 
      xsave=x
182
 
      q2save=qsq
183
 
      if(init.ne.0) goto 10
184
 
c        open(unit=33,file='mrst2002nnlo.dat',status='old')
185
 
        call OpenData(Tablefile)
186
 
        do 20 n=1,nx-1
187
 
        do 20 m=1,nq
188
 
        read(IU,50)f1(n,m),f2(n,m),f3(n,m),f4(n,m),
189
 
     .            f5(n,m),f7(n,m),f6(n,m),f8(n,m)
190
 
c notation: 1=uval 2=val 3=glue 4=usea 5=chm 6=str 7=btm 8=dsea
191
 
  20  continue
192
 
      do 40 m=1,nq
193
 
      f1(nx,m)=0.d0
194
 
      f2(nx,m)=0.d0
195
 
      f3(nx,m)=0.d0
196
 
      f4(nx,m)=0.d0
197
 
      f5(nx,m)=0.d0
198
 
      f6(nx,m)=0.d0
199
 
      f7(nx,m)=0.d0
200
 
      f8(nx,m)=0.d0
201
 
  40  continue
202
 
      do n=1,nx
203
 
      xxl(n)=dlog(xx(n))
204
 
      enddo
205
 
      do m=1,nq
206
 
      qql(m)=dlog(qq(m))
207
 
      enddo
208
 
 
209
 
      call jeppe1(nx,nq,xxl,qql,f1,cc1)
210
 
      call jeppe1(nx,nq,xxl,qql,f2,cc2)
211
 
      call jeppe1(nx,nq,xxl,qql,f3,cc3)
212
 
      call jeppe1(nx,nq,xxl,qql,f4,cc4)
213
 
      call jeppe1(nx,nq,xxl,qql,f6,cc6)
214
 
      call jeppe1(nx,nq,xxl,qql,f8,cc8)
215
 
 
216
 
      emc2=2.045
217
 
      emb2=18.5
218
 
 
219
 
      do 44 m=1,nqc
220
 
      qqlc(m)=qql(m+nqc0)
221
 
      do 44 n=1,nx
222
 
      fc(n,m)=f5(n,m+nqc0)
223
 
   44 continue
224
 
      qqlc(1)=dlog(emc2)
225
 
      call jeppe1(nx,nqc,xxl,qqlc,fc,ccc)
226
 
 
227
 
      do 45 m=1,nqb
228
 
      qqlb(m)=qql(m+nqb0)
229
 
      do 45 n=1,nx
230
 
      fb(n,m)=f7(n,m+nqb0)
231
 
   45 continue
232
 
      qqlb(1)=dlog(emb2)
233
 
      call jeppe1(nx,nqb,xxl,qqlb,fb,ccb)
234
 
 
235
 
 
236
 
      init=1
237
 
   10 continue
238
 
      
239
 
      xlog=dlog(x)
240
 
      qsqlog=dlog(qsq)
241
 
 
242
 
      call jeppe2(xlog,qsqlog,nx,nq,xxl,qql,cc1,upv)
243
 
      call jeppe2(xlog,qsqlog,nx,nq,xxl,qql,cc2,dnv)
244
 
      call jeppe2(xlog,qsqlog,nx,nq,xxl,qql,cc3,glu)
245
 
      call jeppe2(xlog,qsqlog,nx,nq,xxl,qql,cc4,usea)
246
 
      call jeppe2(xlog,qsqlog,nx,nq,xxl,qql,cc6,str)
247
 
      call jeppe2(xlog,qsqlog,nx,nq,xxl,qql,cc8,dsea)
248
 
 
249
 
      chm=0.d0
250
 
      if(qsq.gt.emc2) then 
251
 
      call jeppe2(xlog,qsqlog,nx,nqc,xxl,qqlc,ccc,chm)
252
 
      endif
253
 
 
254
 
      bot=0.d0
255
 
      if(qsq.gt.emb2) then 
256
 
      call jeppe2(xlog,qsqlog,nx,nqb,xxl,qqlb,ccb,bot)
257
 
      endif
258
 
 
259
 
      x=xsave
260
 
      qsq=q2save
261
 
      return
262
 
   50 format(8f10.5)
263
 
      end
264