← Back to branch summary

~maddevelopers/mg5amcnlo/WWW5_caching

« back to all changes in this revision

Viewing changes to users/mardelcourt/PROC_141512/PROC_141512/Source/PDF/pdf.f

  • Committer: John Doe
  • Date: 2013-03-25 20:27:02 UTC
  • Revision ID: john.doe@gmail.com-20130325202702-5sk3t1r8h33ca4p4
first clean version

Show diffs side-by-side

added added

removed removed

Lines of Context:
users/mardelcourt/PROC_141512/PROC_141512/Source/PDF/pdf.f
 
1
      subroutine pftopdg(ih,x,q,pdf)
 
2
c***************************************************************************
 
3
c     Wrapper for calling the pdf of MCFM
 
4
c***************************************************************************
 
5
      implicit none
 
6
c
 
7
c     Arguments
 
8
c
 
9
      DOUBLE  PRECISION x,q,pdf(-7:7)
 
10
      INTEGER IH
 
11
C
 
12
C     Include
 
13
C
 
14
      include 'pdf.inc'
 
15
C      
 
16
      call fdist(ih,x, q, pdf)
 
17
      
 
18
      return    
 
19
      end
 
20
 
 
21
      subroutine fdist(ih,x,xmu,fx)
 
22
C***********************************************************************
 
23
C     MCFM PDF CALLING ROUTINE
 
24
C***********************************************************************
 
25
      implicit none
 
26
      integer ih
 
27
      double precision fx(-7:7),x,xmu
 
28
      double precision u_val,d_val,u_sea,d_sea,s_sea,c_sea,b_sea,gluon
 
29
      double precision Ctq3df,Ctq4Fn,Ctq5Pdf,Ctq6Pdf,Ctq5L
 
30
      double precision q2max
 
31
      double precision epa_electron,epa_proton
 
32
      include 'pdf.inc'
 
33
 
 
34
      integer mode,Iprtn,Irt
 
35
 
 
36
          do Iprtn=-7,7
 
37
             fx(Iprtn)=0d0
 
38
          enddo
 
39
C---set to zero if x out of range
 
40
      if (x .ge. 1d0) then
 
41
          return
 
42
      endif
 
43
                 
 
44
      if     ((pdlabel(1:3) .eq. 'mrs')
 
45
     .   .or. (pdlabel(2:4) .eq. 'mrs')) then
 
46
 
 
47
             if     (pdlabel .eq. 'mrs02nl') then
 
48
             mode=1
 
49
             call mrst2002(x,xmu,mode,u_val,d_val,u_sea,d_sea,
 
50
     &                          s_sea,c_sea,b_sea,gluon)
 
51
             elseif     (pdlabel .eq. 'mrs02nn') then
 
52
             mode=2
 
53
             call mrst2002(x,xmu,mode,u_val,d_val,u_sea,d_sea,
 
54
     &                          s_sea,c_sea,b_sea,gluon)
 
55
             elseif     (pdlabel .eq. 'mrs0119') then
 
56
             mode=1
 
57
             call mrst2001(x,xmu,mode,u_val,d_val,u_sea,d_sea,
 
58
     &                          s_sea,c_sea,b_sea,gluon)
 
59
             elseif (pdlabel .eq. 'mrs0117') then
 
60
             mode=2
 
61
             call mrst2001(x,xmu,mode,u_val,d_val,u_sea,d_sea,
 
62
     &                          s_sea,c_sea,b_sea,gluon)
 
63
             elseif (pdlabel .eq. 'mrs0121') then
 
64
             mode=3
 
65
             call mrst2001(x,xmu,mode,u_val,d_val,u_sea,d_sea,
 
66
     &                          s_sea,c_sea,b_sea,gluon)
 
67
             elseif (pdlabel .eq. 'mrs01_j') then
 
68
             mode=4
 
69
             call mrst2001(x,xmu,mode,u_val,d_val,u_sea,d_sea,
 
70
     &                          s_sea,c_sea,b_sea,gluon)
 
71
             elseif     (pdlabel .eq. 'mrs99_1') then
 
72
             mode=1
 
73
             call mrs99(x,xmu,mode,u_val,d_val,u_sea,d_sea,
 
74
     &                          s_sea,c_sea,b_sea,gluon)
 
75
             elseif (pdlabel .eq. 'mrs99_2') then
 
76
             mode=2
 
77
             call mrs99(x,xmu,mode,u_val,d_val,u_sea,d_sea,
 
78
     &                          s_sea,c_sea,b_sea,gluon)
 
79
             elseif (pdlabel .eq. 'mrs99_3') then
 
80
             mode=3
 
81
             call mrs99(x,xmu,mode,u_val,d_val,u_sea,d_sea,
 
82
     &                          s_sea,c_sea,b_sea,gluon)
 
83
             elseif (pdlabel .eq. 'mrs99_4') then
 
84
             mode=4
 
85
             call mrs99(x,xmu,mode,u_val,d_val,u_sea,d_sea,
 
86
     &                          s_sea,c_sea,b_sea,gluon)
 
87
             elseif (pdlabel .eq. 'mrs99_5') then
 
88
             mode=5
 
89
             call mrs99(x,xmu,mode,u_val,d_val,u_sea,d_sea,
 
90
     &                          s_sea,c_sea,b_sea,gluon)
 
91
             elseif (pdlabel .eq. 'mrs99_6') then
 
92
             mode=6
 
93
             call mrs99(x,xmu,mode,u_val,d_val,u_sea,d_sea,
 
94
     &                          s_sea,c_sea,b_sea,gluon)
 
95
             elseif (pdlabel .eq. 'mrs99_7') then
 
96
             mode=7
 
97
             call mrs99(x,xmu,mode,u_val,d_val,u_sea,d_sea,
 
98
     &                          s_sea,c_sea,b_sea,gluon)
 
99
             elseif (pdlabel .eq. 'mrs99_8') then
 
100
             mode=8
 
101
             call mrs99(x,xmu,mode,u_val,d_val,u_sea,d_sea,
 
102
     &                          s_sea,c_sea,b_sea,gluon)
 
103
             elseif (pdlabel .eq. 'mrs99_9') then
 
104
             mode=9
 
105
             call mrs99(x,xmu,mode,u_val,d_val,u_sea,d_sea,
 
106
     &                          s_sea,c_sea,b_sea,gluon)
 
107
             elseif (pdlabel .eq. 'mrs9910') then
 
108
             mode=10
 
109
             call mrs99(x,xmu,mode,u_val,d_val,u_sea,d_sea,
 
110
     &                          s_sea,c_sea,b_sea,gluon)
 
111
             elseif (pdlabel .eq. 'mrs9911') then
 
112
             mode=11
 
113
             call mrs99(x,xmu,mode,u_val,d_val,u_sea,d_sea,
 
114
     &                          s_sea,c_sea,b_sea,gluon)
 
115
             elseif (pdlabel .eq. 'mrs9912') then
 
116
             mode=12
 
117
             call mrs99(x,xmu,mode,u_val,d_val,u_sea,d_sea,
 
118
     &                          s_sea,c_sea,b_sea,gluon)
 
119
             elseif (pdlabel .eq. 'mrs98z1') then
 
120
             mode=1
 
121
             call mrs98(x,xmu,mode,u_val,d_val,u_sea,d_sea,
 
122
     &                          s_sea,c_sea,b_sea,gluon)
 
123
             elseif (pdlabel .eq. 'mrs98z2') then
 
124
             mode=2 
 
125
             call mrs98(x,xmu,mode,u_val,d_val,u_sea,d_sea,
 
126
     &                          s_sea,c_sea,b_sea,gluon)
 
127
             elseif (pdlabel .eq. 'mrs98z3') then
 
128
             mode=3
 
129
             call mrs98(x,xmu,mode,u_val,d_val,u_sea,d_sea,
 
130
     &                          s_sea,c_sea,b_sea,gluon)
 
131
             elseif (pdlabel .eq. 'mrs98z4') then
 
132
             mode=4
 
133
             call mrs98(x,xmu,mode,u_val,d_val,u_sea,d_sea,
 
134
     &                          s_sea,c_sea,b_sea,gluon)
 
135
             elseif (pdlabel .eq. 'mrs98z5') then
 
136
             mode=5
 
137
             call mrs98(x,xmu,mode,u_val,d_val,u_sea,d_sea,
 
138
     &                          s_sea,c_sea,b_sea,gluon)
 
139
             elseif (pdlabel .eq. 'mrs98l1') then
 
140
             mode=1
 
141
             call mrs98lo(x,xmu,mode,u_val,d_val,u_sea,d_sea,
 
142
     &                          s_sea,c_sea,b_sea,gluon)
 
143
             elseif (pdlabel .eq. 'mrs98l2') then
 
144
             mode=2 
 
145
             call mrs98lo(x,xmu,mode,u_val,d_val,u_sea,d_sea,
 
146
     &                          s_sea,c_sea,b_sea,gluon)
 
147
             elseif (pdlabel .eq. 'mrs98l3') then
 
148
             mode=3
 
149
             call mrs98lo(x,xmu,mode,u_val,d_val,u_sea,d_sea,
 
150
     &                          s_sea,c_sea,b_sea,gluon)
 
151
             elseif (pdlabel .eq. 'mrs98l4') then
 
152
             mode=4
 
153
             call mrs98lo(x,xmu,mode,u_val,d_val,u_sea,d_sea,
 
154
     &                          s_sea,c_sea,b_sea,gluon)
 
155
             elseif (pdlabel .eq. 'mrs98l5') then
 
156
             mode=5
 
157
             call mrs98lo(x,xmu,mode,u_val,d_val,u_sea,d_sea,
 
158
     &                          s_sea,c_sea,b_sea,gluon)
 
159
             elseif (pdlabel .eq. 'mrs98ht') then
 
160
             mode=1
 
161
             call mrs98ht(x,xmu,mode,u_val,d_val,u_sea,d_sea,
 
162
     &                          s_sea,c_sea,b_sea,gluon)
 
163
             endif
 
164
c-----assign mrs to standard grid
 
165
            fx(-5)=b_sea/x
 
166
            fx(-4)=c_sea/x
 
167
            fx(-3)=s_sea/x
 
168
            fx( 0)=gluon/x
 
169
            fx(+3)=fx(-3)
 
170
            fx(+4)=fx(-4)
 
171
            fx(+5)=fx(-5)
 
172
               fx(1)=(d_val+d_sea)/x
 
173
               fx(2)=(u_val+u_sea)/x
 
174
               fx(-1)=d_sea/x
 
175
               fx(-2)=u_sea/x
 
176
C
 
177
      elseif (pdlabel(1:5) .eq. 'cteq3') then
 
178
C     
 
179
         if (pdlabel .eq. 'cteq3_m') then
 
180
            mode=1
 
181
         elseif (pdlabel .eq. 'cteq3_l') then
 
182
            mode=2
 
183
         elseif (pdlabel .eq. 'cteq3_d') then
 
184
            mode=3
 
185
         endif
 
186
         fx(-5)=Ctq3df(mode,-5,x,xmu,Irt)/x
 
187
         fx(-4)=Ctq3df(mode,-4,x,xmu,Irt)/x
 
188
         fx(-3)=Ctq3df(mode,-3,x,xmu,Irt)/x
 
189
         
 
190
         fx(0)=Ctq3df(mode,0,x,xmu,Irt)/x
 
191
         
 
192
         fx(+3)=Ctq3df(mode,+3,x,xmu,Irt)/x
 
193
         fx(+4)=Ctq3df(mode,+4,x,xmu,Irt)/x
 
194
         fx(+5)=Ctq3df(mode,+5,x,xmu,Irt)/x
 
195
            fx(-1)=Ctq3df(mode,-2,x,xmu,Irt)/x
 
196
            fx(-2)=Ctq3df(mode,-1,x,xmu,Irt)/x
 
197
            fx(1)=Ctq3df(mode,+2,x,xmu,Irt)/x+fx(-1)
 
198
            fx(2)=Ctq3df(mode,+1,x,xmu,Irt)/x+fx(-2)
 
199
C     
 
200
      elseif (pdlabel(1:5) .eq. 'cteq4') then
 
201
C     
 
202
         if (pdlabel .eq. 'cteq4_m') then
 
203
            mode=1
 
204
         elseif (pdlabel .eq. 'cteq4_d') then
 
205
            mode=2
 
206
         elseif (pdlabel .eq. 'cteq4_l') then
 
207
            mode=3
 
208
         elseif (pdlabel .eq. 'cteq4a1') then
 
209
            mode=4
 
210
         elseif (pdlabel .eq. 'cteq4a2') then
 
211
            mode=5
 
212
         elseif (pdlabel .eq. 'cteq4a3') then
 
213
            mode=6
 
214
         elseif (pdlabel .eq. 'cteq4a4') then
 
215
            mode=7
 
216
         elseif (pdlabel .eq. 'cteq4a5') then
 
217
            mode=8
 
218
         elseif (pdlabel .eq. 'cteq4hj') then
 
219
            mode=9
 
220
         elseif (pdlabel .eq. 'cteq4lq') then
 
221
            mode=10
 
222
         endif
 
223
         
 
224
         fx(-5)=Ctq4Fn(mode,-5,x,xmu)
 
225
         fx(-4)=Ctq4Fn(mode,-4,x,xmu)
 
226
         fx(-3)=Ctq4Fn(mode,-3,x,xmu)
 
227
         
 
228
         fx(0)=Ctq4Fn(mode,0,x,xmu)
 
229
         
 
230
         fx(+3)=Ctq4Fn(mode,+3,x,xmu)
 
231
         fx(+4)=Ctq4Fn(mode,+4,x,xmu)
 
232
         fx(+5)=Ctq4Fn(mode,+5,x,xmu)
 
233
            fx(1)=Ctq4Fn(mode,+2,x,xmu)
 
234
            fx(2)=Ctq4Fn(mode,+1,x,xmu)
 
235
            fx(-1)=Ctq4Fn(mode,-2,x,xmu)
 
236
            fx(-2)=Ctq4Fn(mode,-1,x,xmu)
 
237
C
 
238
      elseif (pdlabel .eq. 'cteq5l1') then
 
239
C
 
240
         fx(-5)=Ctq5L(-5,x,xmu)
 
241
         fx(-4)=Ctq5L(-4,x,xmu)
 
242
         fx(-3)=Ctq5L(-3,x,xmu)
 
243
         
 
244
         fx(0)=Ctq5L(0,x,xmu)
 
245
         
 
246
         fx(+3)=Ctq5L(+3,x,xmu)
 
247
         fx(+4)=Ctq5L(+4,x,xmu)
 
248
         fx(+5)=Ctq5L(+5,x,xmu)
 
249
         
 
250
            fx(1)=Ctq5L(+2,x,xmu)
 
251
            fx(2)=Ctq5L(+1,x,xmu)
 
252
            fx(-1)=Ctq5L(-2,x,xmu)
 
253
            fx(-2)=Ctq5L(-1,x,xmu)
 
254
C         
 
255
      elseif ((pdlabel(1:5) .eq. 'cteq5') .or. 
 
256
     .        (pdlabel(1:4) .eq. 'ctq5')) then
 
257
C         
 
258
         fx(-5)=Ctq5Pdf(-5,x,xmu)
 
259
         fx(-4)=Ctq5Pdf(-4,x,xmu)
 
260
         fx(-3)=Ctq5Pdf(-3,x,xmu)
 
261
         
 
262
         fx(0)=Ctq5Pdf(0,x,xmu)
 
263
         
 
264
         fx(+3)=Ctq5Pdf(+3,x,xmu)
 
265
         fx(+4)=Ctq5Pdf(+4,x,xmu)
 
266
         fx(+5)=Ctq5Pdf(+5,x,xmu)
 
267
         
 
268
            fx(1)=Ctq5Pdf(+2,x,xmu)
 
269
            fx(2)=Ctq5Pdf(+1,x,xmu)
 
270
            fx(-1)=Ctq5Pdf(-2,x,xmu)
 
271
            fx(-2)=Ctq5Pdf(-1,x,xmu)
 
272
C                  
 
273
      elseif (pdlabel(1:5) .eq. 'cteq6') then
 
274
C         
 
275
         fx(-5)=Ctq6Pdf(-5,x,xmu)
 
276
         fx(-4)=Ctq6Pdf(-4,x,xmu)
 
277
         fx(-3)=Ctq6Pdf(-3,x,xmu)
 
278
         
 
279
         fx(0)=Ctq6Pdf(0,x,xmu)
 
280
         
 
281
         fx(+3)=Ctq6Pdf(+3,x,xmu)
 
282
         fx(+4)=Ctq6Pdf(+4,x,xmu)
 
283
         fx(+5)=Ctq6Pdf(+5,x,xmu)
 
284
         
 
285
            fx(1)=Ctq6Pdf(+2,x,xmu)
 
286
            fx(2)=Ctq6Pdf(+1,x,xmu)
 
287
            fx(-1)=Ctq6Pdf(-2,x,xmu)
 
288
            fx(-2)=Ctq6Pdf(-1,x,xmu)
 
289
      endif      
 
290
c
 
291
c  a "diffractive" photon
 
292
c      
 
293
      q2max=xmu*xmu
 
294
      if(ih .eq. 3) then  !from the electron
 
295
          fx(7)=epa_electron(x,q2max)
 
296
      elseif(ih .eq. 2) then  !from a proton without breaking
 
297
          fx(7)=epa_proton(x,q2max)
 
298
      endif      
 
299
      
 
300
      return
 
301
      end
 
302
      
 
303
  
 
304