]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PYTHIA6/lhaglue.F
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
Constructors corrected.
[u/mrichter/AliRoot.git] / PYTHIA6 / lhaglue.F
CommitLineData
48f2b607 1C*********************************************************************
2
3C... LHAGLUE Interface to LHAPDF library of modern parton
4C... density functions (PDF) with uncertainties
5C...
6C...Authors for v4: Dimitri Bourilkov, Craig Group, Mike Whalley
7C...
8C...Authors for v3: Dimitri Bourilkov, Craig Group, Mike Whalley
9C...
10C...Author for v1 and v2: Dimitri Bourilkov bourilkov@mailaps.org
11C... University of Florida
12C...
13C...HERWIG interface by Dimitri Bourilkov and Craig Group
14C...
15C...New numbering scheme and upgrade for LHAPDF v2.1
16C...by Dimitri Bourilkov and Mike Whalley
17C...
18C...For more information, or when you cite this interface, currently
19C...the official reference is:
20C...D.Bourilkov, "Study of Parton Density Function Uncertainties with
21C...LHAPDF and PYTHIA at LHC", hep-ph/0305126.
22C...
23C...The official LHAPDF page is:
24C...
25C... http://durpdg.dur.ac.uk/lhapdf/index.html
26C...
27C...The interface contains four subroutines (similar to PDFLIB).
28C...It can be used seamlessly by Monte Carlo generators
29C...interfaced to PDFLIB or in stand-alone mode.
30C...
31C... For initialization (called once)
32C...
33C... PDFSET(PARM,VALUE)
34C...
35C... For the proton/pion structure functions
36C...
37C... STRUCTM(X,Q,UPV,DNV,USEA,DSEA,STR,CHM,BOT,TOP,GLU)
38C...
39C... For the photon structure functions
40C...
41C... STRUCTP(X,Q2,P2,IP2,UPV,DNV,USEA,DSEA,STR,CHM,BOT,TOP,GLU)
42C...
43C... For statistics ON structure functions (under/over-flows)
44C...
45C... PDFSTA
46C...
47C...This interface can be invoked in 3 ways depending
48C...on the value of PARM(1) provided by the user when
49C...calling PDFSET(PARM,VALUE):
50C...
51C... For PYTHIA: PARM(1).EQ.'NPTYPE'
52C... (this is set automatically by PYTHIA)
53C...
54C... For HERWIG: PARM(1).EQ.'HWLHAPDF'
55C... (set by the USER e.g. in the main program like this:
56C... AUTPDF(1) = 'HWLHAPDF'
57C... AUTPDF(2) = 'HWLHAPDF' )
58C...
59C... For Stand-alone: PARM(1).EQ.'DEFAULT'
60C... (can be used for PDF studies or when interfacing
61C... new generators)
62C...
63C...The LHAPDF set/member is selected depending on the value of:
64C...
65C... PYTHIA: ABS(MSTP(51)) - proton
66C... ABS(MSTP(53)) - pion
67C... ABS(MSTP(55)) - photon
68C...
69C... HERWIG: ABS(INT(VALUE(1)))
70C...
71C... STAND-ALONE: ABS(INT(VALUE(1)))
72C...
73C...
74C... CONTROL switches:
75C... ==================
76C...
77C... THE LOCATION OF THE LHAPDF LIBRARY HAS TO BE SPECIFIED
78C... AS DESCRIBED BELOW (the rest is optional)
79C...
80C... if the user does nothing, sensible defaults
81C... are active; to change the behaviour, the corresponding
82C... values of LHAPARM() should be set to the values given below
83C...
84C... Location of the LHAPDF library of PDFs (pathname):
85C... uses common block /LHAPDFC/
86C...
87C... If the user does nothing => default = subdir PDFsets of the
88C... current directory (can be real subdir
89C... OR a soft link to the real location)
90C... If the user sets LHAPATH => supplied by the USER who defines the
91C... path in common block COMMON/LHAPDFC/LHAPATH
92C... BEFORE calling PDFSET
93C...
94C... Other controls:
95C... ===============
96C... use common block /LHACONTROL/
97C...
98C... Collect statistics on under/over-flow requests for PDFs
99C... outside their validity ranges in X and Q**2
100C... (call PDFSTA at end of run to print it out)
101C...
102C... LHAPARM(16).EQ.'NOSTAT' => No statistics (faster)
103C... LHAPARM(16).NE.'NOSTAT' => Default: collect statistics
104C...
105C... Option to use the values for the strong coupling alpha_s
106C... as computed in LHAPDF in the MC generator
107C... (to ensure uniformity between the MC generator and the PDF set)
108C... WARNING: implemented ONLY for PYTHIA in LHAPDFv4
109C...
110C... LHAPARM(17).EQ.'LHAPDF' => Use alpha_s from LHAPDF
111C... LHAPARM(17).NE.'LHAPDF' => Default (same as LHAPDF v1/v3)
112C...
113C... Extrapolation of PDFs outside LHAPDF validity range given by
114C... [Xmin,Xmax] and [Q2min,Q2max]; DEFAULT => PDFs "frozen" at the
115C... boundaries
116C...
117C... LHAPARM(18).EQ.'EXTRAPOLATE' => Extrapolate PDFs on OWN RISK
118C... WARNING: Crazy values can be returned
119C...
120C... Printout of initialization information in PDFSET (by default)
121C...
122C... LHAPARM(19).EQ.'SILENT' => No printout (silent mode)
123C... LHAPARM(19).EQ.'LOWKEY' => Print 5 times (almost silent mode)
124C...
125C...
126C...v5.0 06-Oct-2005 Major change to allow multiset-initializations
127C...v4.0 28-Apr-2005 PDFSTA routine; option to use Alfa_s from LHAPDF
128C...v4.0 21-Mar-2005 Photon/pion/new p PDFs, updated for LHAPDF v4
129C...v3.1 26-Apr-2004 New numbering scheme, updated for LHAPDF v2/v3
130C...v3.0 23-Jan-2004 HERWIG interface added
131C...v2.0 20-Sep-2003 PDFLIB style adopted
132C...v1.0 05-Mar-2003 First working version from PYTHIA to LHAPDF v1
133C...
134C...interface to LHAPDF library
135
136C*********************************************************************
137
138C...PDFSET
139C...Initialization for use of parton distributions
140C... according to the LHAPDF interface.
141C...
142C...v4.0 28-Apr-2005 Option to use Alfa_s from LHAPDF
143C...v4.0 21-Mar-2005 Photon/pion/new p PDFs, updated for LHAPDF v4
144C...v3.1 26-Apr-2004 New numbering scheme
145C...v3.0 23-Jan-2004 HERWIG interface added
146C...
147C...interface to LHAPDF library
148
149 SUBROUTINE PDFSET(PARM,VALUE)
150C...Double precision and integer declarations.
151 IMPLICIT DOUBLE PRECISION(A-H, O-Z)
152 IMPLICIT INTEGER(I-N)
153c...additions for multiset use
154 integer noemax,nopmax,npfmax,nofmax,linemax,nmxset
155 parameter (noemax=1000,nopmax=40,npfmax=10,nofmax=10,linemax=20)
156c nmxset is the max number of sets that can be initialised at one time ---- added V5
157 parameter (nmxset=3)
158
159 character*232 lhapath
160 character*272 lhaname,lhanames(nmxset)
161c character*172 LHANAMES(nmxset)
162 integer LHAMEMBERS(nmxset),LHANUMBERS(nmxset)
163 common/LHASETS/LHANAMES,LHANUMBERS,LHAMEMBERS,nsets
164 real*8 xxmin(nmxset),xxmax(nmxset),qq2min(nmxset),qq2max(nmxset)
165 save xxmin,xxmax,qq2min,qq2max
166C...Commonblocks.
167 COMMON/PYDAT1/MSTU(200),PARU(200),MSTJ(200),PARJ(200)
168 SAVE /PYDAT1/
169 COMMON/PYPARS/MSTP(200),PARP(200),MSTI(200),PARI(200)
170 SAVE /PYPARS/
171c... following 2 for earlier Pythia versions
172 COMMON/LUDAT1/MSTU5(200),PARU5(200),MSTJ5(200),PARJ5(200)
173 SAVE /LUDAT1/
174C...Interface to LHAPDFLIB.
175c CHARACTER*172 LHANAME
176 INTEGER LHASET, LHAMEMB
177 COMMON/LHAPDF/LHANAME, LHASET, LHAMEMB
178 SAVE /LHAPDF/
179 DOUBLE PRECISION QCDLHA4, QCDLHA5
180 INTEGER NFLLHA
181 COMMON/LHAPDFR/QCDLHA4, QCDLHA5, NFLLHA
182 SAVE /LHAPDFR/
183c CHARACTER*132 LHAPATH
184 COMMON/LHAPDFC/LHAPATH
185 SAVE /LHAPDFC/
186 CHARACTER*20 LHAPARM(20)
187 DOUBLE PRECISION LHAVALUE(20)
188 COMMON/LHACONTROL/LHAPARM,LHAVALUE
189 SAVE/LHACONTROL/
190 INTEGER LHAEXTRP
191 COMMON/LHAPDFE/LHAEXTRP
192 SAVE /LHAPDFE/
193 INTEGER LHASILENT
194 COMMON/LHASILENT/LHASILENT
195 SAVE /LHASILENT/
196 DOUBLE PRECISION XMINNUM,XMAXNUM,Q2MINNUM,Q2MAXNUM,TOTNUM,
197 > XMINNUP,XMAXNUP,Q2MINNUP,Q2MAXNUP,TOTNUP
198 COMMON/LHAGLSTA/ XMINNUM,XMAXNUM,Q2MINNUM,Q2MAXNUM,TOTNUM,
199 > XMINNUP,XMAXNUP,Q2MINNUP,Q2MAXNUP,TOTNUP
200 SAVE/LHAGLSTA/
201C...Interface to PDFLIB.
202 COMMON/W50511/ NPTYPEPDFL,NGROUPPDFL,NSETPDFL,MODEPDFL,
203 > NFLPDFL,LOPDFL,TMASPDFL
204 SAVE /W50511/
205 DOUBLE PRECISION TMASPDFL
206C...Interface to PDFLIB.
207 COMMON/W50512/QCDL4,QCDL5
208 SAVE /W50512/
209 DOUBLE PRECISION QCDL4,QCDL5
210C...Interface to PDFLIB.
211 COMMON/W50513/XMIN,XMAX,Q2MIN,Q2MAX
212 SAVE /W50513/
213 DOUBLE PRECISION XMIN,XMAX,Q2MIN,Q2MAX
214C...Local arrays and character variables (NOT USED here DB)
215 CHARACTER*20 PARM(20)
216 DOUBLE PRECISION VALUE(20)
217 INTEGER LHAPATHLEN
218 INTEGER LHAINPUT
219 INTEGER LHASELECT
220 INTEGER LHAPRINT
221 INTEGER LHAONCE
222 INTEGER LHAFIVE
223 SAVE LHAONCE
224 SAVE LHAFIVE
225 logical first
226 save first
227
228 INTEGER LNROOT
229 CHARACTER*1000 CHROOT
230 CHROOT=' '
231 DATA LHAONCE/0/
232 DATA LHAFIVE/0/
233 data first/.TRUE./
234
235 if(first .AND. (LHAPARM(20).NE.'LHAPATH')) then
236c...overide the default PDFsets path
237c ... check first if the environmental variable LHAPATH is set
238 call getenv('LHAPATH',lhapath)
239 if(lhapath.eq.'') then
240C The environment variable LHAPATH is not set.
241C Take the data from $ALICE_ROOT/LHAPDF/PDFsets
242 CALL GETENV('ALICE_ROOT',CHROOT)
243 LNROOT = LNBLNK(CHROOT)
244 IF(LNROOT.LE.0) THEN
245 LHAPATH='PDFsets' ! Default value
246 ELSE
247 LHAPATH=CHROOT(1:LNROOT)//'/LHAPDF/PDFsets'
248 ENDIF
249 endif
250 first=.FALSE.
251 endif
252c
253*
254C...Init
255 LHAEXTRP = 0
256 IF(LHAPARM(18).EQ.'EXTRAPOLATE')
257 > THEN ! Extrapolate PDFs on own risk
258 LHAEXTRP = 1
259 ENDIF
260 LHASILENT = 0
261 IF(LHAPARM(19).EQ.'SILENT') THEN ! No printout (silent MODE)
262 LHASILENT = 1
263 ELSEIF(LHAPARM(19).EQ.'LOWKEY') THEN ! Print 5 times (lowkey MODE)
264 IF(LHAFIVE .LT. 6) THEN
265 LHAFIVE = LHAFIVE + 1
266 ELSE
267 LHASILENT = 1
268 ENDIF
269 ENDIF
270 IF(PARM(1).EQ.'NPTYPE') THEN ! PYTHIA
271 if(MSTP(181).ge.6) then
272 LHAPRINT = MSTU(11)
273 else
274 LHAPRINT = MSTU5(11)
275 endif
276 IF(VALUE(1) .EQ. 1) THEN ! nucleon
277 LHAINPUT = ABS(MSTP(51))
278 ELSEIF(VALUE(1) .EQ. 2) THEN ! pion
279 LHAINPUT = ABS(MSTP(53))
280 ELSEIF(VALUE(1) .EQ. 3) THEN ! photon
281 LHAINPUT = ABS(MSTP(55))
282 ENDIF
283 IF(LHASILENT .NE. 1)
284 > PRINT *,'==== PYTHIA WILL USE LHAPDF ===='
285 ELSEIF(PARM(1).EQ.'HWLHAPDF') THEN ! HERWIG
286 LHAINPUT = ABS(INT(VALUE(1)))
287 IF(LHAONCE.EQ.LHAINPUT) RETURN
288 IF(LHASILENT .NE. 1)
289 > PRINT *,'==== HERWIG WILL USE LHAPDF ===='
290 LHAPRINT = 6
291 LHAONCE = LHAINPUT
292 ELSEIF(PARM(1).EQ.'DEFAULT') THEN ! Stand-alone
293 LHAINPUT = ABS(INT(VALUE(1)))
294 IF(LHAONCE.EQ.LHAINPUT) RETURN
295 IF(LHASILENT .NE. 1)
296 > PRINT *,'==== STAND-ALONE LHAGLUE MODE TO USE LHAPDF ===='
297 LHAPRINT = 6
298 LHAONCE = LHAINPUT
299 ELSE
300 PRINT *,'== UNKNOWN LHAPDF INTERFACE CALL! STOP EXECUTION! =='
301 STOP
302 ENDIF
303C...Initialize parton distributions: LHAPDFLIB.
304 LHAPATHLEN=INDEX(LHAPATH,' ')-1
305 LHASET = LHAINPUT
306 XMIN = 1.0D-6 ! X_min for current PDF set
307 XMAX = 1.0D0 ! X_max for current PDF set
308 Q2MIN = 1.0D0**2 ! Q**2_min scale for current PDF set [GeV]
309 Q2MAX = 1.0D5**2 ! Q**2_max scale for current PDF set [GeV]
310C...
311C...Protons
312C...
313C...CTEQ Family
314 IF((LHAINPUT .GE. 10000) .AND. (LHAINPUT .LE. 19999)) THEN
315 Q2MAX = 1.0D08
316 IF((LHAINPUT .GE. 10000) .AND. (LHAINPUT .LE. 10040)) THEN
317 LHASET = 10000
318 LHANAME=LHAPATH(1:LHAPATHLEN)//'/cteq6.LHpdf'
319 Q2MIN = 1.69D0
320 ELSEIF((LHAINPUT .GE. 10041) .AND. (LHAINPUT .LE. 10041)) THEN
321 LHASET = 10041
322 LHANAME=LHAPATH(1:LHAPATHLEN)//'/cteq6l.LHpdf'
323 Q2MIN = 1.69D0
324 ELSEIF((LHAINPUT .GE. 10042) .AND. (LHAINPUT .LE. 10042)) THEN
325 LHASET = 10042
326 LHANAME=LHAPATH(1:LHAPATHLEN)//'/cteq6ll.LHpdf'
327 Q2MIN = 1.69D0
328 ELSEIF((LHAINPUT .GE. 10050) .AND. (LHAINPUT .LE. 10090)) THEN
329 LHASET = 10050
330 LHANAME=LHAPATH(1:LHAPATHLEN)//'/cteq6mE.LHgrid'
331 Q2MIN = 1.69D0
332 ELSEIF((LHAINPUT .GE. 10100) .AND. (LHAINPUT .LE. 10140)) THEN
333 LHASET = 10100
334 LHANAME=LHAPATH(1:LHAPATHLEN)//'/cteq61.LHpdf'
335 Q2MIN = 1.69D0
336 ELSEIF((LHAINPUT .GE. 10150) .AND. (LHAINPUT .LE. 10190)) THEN
337 LHASET = 10150
338 LHANAME=LHAPATH(1:LHAPATHLEN)//'/cteq61.LHgrid'
339 Q2MIN = 1.69D0
340 ELSEIF((LHAINPUT .GE. 10250) .AND. (LHAINPUT .LE. 10269)) THEN
341 LHASET = 10250
342 LHANAME=LHAPATH(1:LHAPATHLEN)//'/cteq6AB.LHgrid'
343 Q2MIN = 1.69D0
344 ELSEIF((LHAINPUT .GE. 19050) .AND. (LHAINPUT .LE. 19050)) THEN
345 LHASET = 19050
346 LHANAME=LHAPATH(1:LHAPATHLEN)//'/cteq5m.LHgrid'
347 XMIN=1.0D-5
348 ELSEIF((LHAINPUT .GE. 19051) .AND. (LHAINPUT .LE. 19051)) THEN
349 LHASET = 19051
350 LHANAME=LHAPATH(1:LHAPATHLEN)//'/cteq5m1.LHgrid'
351 XMIN=1.0D-5
352 ELSEIF((LHAINPUT .GE. 19060) .AND. (LHAINPUT .LE. 19060)) THEN
353 LHASET = 19060
354 LHANAME=LHAPATH(1:LHAPATHLEN)//'/cteq5d.LHgrid'
355 XMIN=1.0D-5
356 ELSEIF((LHAINPUT .GE. 19070) .AND. (LHAINPUT .LE. 19070)) THEN
357 LHASET = 19070
358 LHANAME=LHAPATH(1:LHAPATHLEN)//'/cteq5l.LHgrid'
359 XMIN=1.0D-5
360 ELSEIF((LHAINPUT .GE. 19150) .AND. (LHAINPUT .LE. 19150)) THEN
361 LHASET = 19150
362 LHANAME=LHAPATH(1:LHAPATHLEN)//'/cteq4m.LHgrid'
363 Q2MIN = 2.56D0
364 XMIN=1.0D-5
365 ELSEIF((LHAINPUT .GE. 19160) .AND. (LHAINPUT .LE. 19160)) THEN
366 LHASET = 19160
367 LHANAME=LHAPATH(1:LHAPATHLEN)//'/cteq4d.LHgrid'
368 Q2MIN = 2.56D0
369 XMIN=1.0D-5
370 ELSEIF((LHAINPUT .GE. 19170) .AND. (LHAINPUT .LE. 19170)) THEN
371 LHASET = 19170
372 LHANAME=LHAPATH(1:LHAPATHLEN)//'/cteq4l.LHgrid'
373 Q2MIN = 2.56D0
374 XMIN=1.0D-5
375 ELSE
376 WRITE(LHAPRINT,5150) LHASET
377 STOP
378 ENDIF
379C...MRST Family
380 ELSEIF((LHAINPUT .GE. 20000) .AND. (LHAINPUT .LE. 29999)) THEN
381 Q2MIN = 1.25D0
382 Q2MAX = 1.0D07
383 XMIN = 1.0D-5
384 IF((LHAINPUT .GE. 20000) .AND. (LHAINPUT .LE. 20004)) THEN
385 LHASET = 20000
386 LHANAME=LHAPATH(1:LHAPATHLEN)//'/MRST2001nlo.LHpdf'
387 ELSEIF((LHAINPUT .GE. 20050) .AND. (LHAINPUT .LE. 20054)) THEN
388 LHASET = 20050
389 LHANAME=LHAPATH(1:LHAPATHLEN)//'/MRST2001nlo.LHgrid'
390 ELSEIF((LHAINPUT .GE. 20060) .AND. (LHAINPUT .LE. 20061)) THEN
391 LHASET = 20060
392 LHANAME=LHAPATH(1:LHAPATHLEN)//'/MRST2001lo.LHgrid'
393 ELSEIF((LHAINPUT .GE. 20070) .AND. (LHAINPUT .LE. 20074)) THEN
394 LHASET = 20070
395 LHANAME=LHAPATH(1:LHAPATHLEN)//'/MRST2001nnlo.LHgrid'
396 ELSEIF((LHAINPUT .GE. 20100) .AND. (LHAINPUT .LE. 20130)) THEN
397 LHASET = 20100
398 LHANAME=LHAPATH(1:LHAPATHLEN)//'/MRST2001E.LHpdf'
399 ELSEIF((LHAINPUT .GE. 20150) .AND. (LHAINPUT .LE. 20180)) THEN
400 LHASET = 20150
401 LHANAME=LHAPATH(1:LHAPATHLEN)//'/MRST2001E.LHgrid'
402 ELSEIF((LHAINPUT .GE. 20200) .AND. (LHAINPUT .LE. 20201)) THEN
403 LHASET = 20200
404 LHANAME=LHAPATH(1:LHAPATHLEN)//'/MRST2002nlo.LHpdf'
405 ELSEIF((LHAINPUT .GE. 20250) .AND. (LHAINPUT .LE. 20251)) THEN
406 LHASET = 20250
407 LHANAME=LHAPATH(1:LHAPATHLEN)//'/MRST2002nlo.LHgrid'
408 ELSEIF((LHAINPUT .GE. 20270) .AND. (LHAINPUT .LE. 20271)) THEN
409 LHASET = 20270
410 LHANAME=LHAPATH(1:LHAPATHLEN)//'/MRST2002nnlo.LHgrid'
411 ELSEIF((LHAINPUT .GE. 20300) .AND. (LHAINPUT .LE. 20301)) THEN
412 LHASET = 20300
413 LHANAME=LHAPATH(1:LHAPATHLEN)//'/MRST2003cnlo.LHpdf'
414 Q2MIN = 10.0D0
415 XMIN = 1.0D-3
416 ELSEIF((LHAINPUT .GE. 20350) .AND. (LHAINPUT .LE. 20351)) THEN
417 LHASET = 20350
418 LHANAME=LHAPATH(1:LHAPATHLEN)//'/MRST2003cnlo.LHgrid'
419 Q2MIN = 10.0D0
420 XMIN = 1.0D-3
421 ELSEIF((LHAINPUT .GE. 20370) .AND. (LHAINPUT .LE. 20371)) THEN
422 LHASET = 20370
423 LHANAME=LHAPATH(1:LHAPATHLEN)//'/MRST2003cnnlo.LHgrid'
424 Q2MIN = 7.0D0
425 XMIN = 1.0D-3
426 ELSEIF((LHAINPUT .GE. 20400) .AND. (LHAINPUT .LE. 20401)) THEN
427 LHASET = 20400
428 LHANAME=LHAPATH(1:LHAPATHLEN)//'/MRST2004nlo.LHpdf'
429 ELSEIF((LHAINPUT .GE. 20406) .AND. (LHAINPUT .LE. 20407)) THEN
430 LHASET = 20406
431 LHANAME=LHAPATH(1:LHAPATHLEN)//'/MRST2004FF3nlo.LHpdf'
432 ELSEIF((LHAINPUT .GE. 20408) .AND. (LHAINPUT .LE. 20409)) THEN
433 LHASET = 20408
434 LHANAME=LHAPATH(1:LHAPATHLEN)//'/MRST2004FF4nlo.LHpdf'
435 ELSEIF((LHAINPUT .GE. 20450) .AND. (LHAINPUT .LE. 20451)) THEN
436 LHASET = 20450
437 LHANAME=LHAPATH(1:LHAPATHLEN)//'/MRST2004nlo.LHgrid'
438 ELSEIF((LHAINPUT .GE. 20452) .AND. (LHAINPUT .LE. 20453)) THEN
439 LHASET = 20452
440 LHANAME=LHAPATH(1:LHAPATHLEN)//'/MRST2004FF3lo.LHgrid'
441 ELSEIF((LHAINPUT .GE. 20454) .AND. (LHAINPUT .LE. 20455)) THEN
442 LHASET = 20454
443 LHANAME=LHAPATH(1:LHAPATHLEN)//'/MRST2004FF4lo.LHgrid'
444 ELSEIF((LHAINPUT .GE. 20456) .AND. (LHAINPUT .LE. 20457)) THEN
445 LHASET = 20456
446 LHANAME=LHAPATH(1:LHAPATHLEN)//'/MRST2004FF3nlo.LHgrid'
447 ELSEIF((LHAINPUT .GE. 20458) .AND. (LHAINPUT .LE. 20459)) THEN
448 LHASET = 20458
449 LHANAME=LHAPATH(1:LHAPATHLEN)//'/MRST2004FF4nlo.LHgrid'
450 ELSEIF((LHAINPUT .GE. 20470) .AND. (LHAINPUT .LE. 20471)) THEN
451 LHASET = 20470
452 LHANAME=LHAPATH(1:LHAPATHLEN)//'/MRST2004nnlo.LHgrid'
453 ELSEIF((LHAINPUT .GE. 29000) .AND. (LHAINPUT .LE. 29003)) THEN
454 LHASET = 29000
455 LHANAME=LHAPATH(1:LHAPATHLEN)//'/MRST98.LHpdf'
456 ELSEIF((LHAINPUT .GE. 29040) .AND. (LHAINPUT .LE. 29045)) THEN
457 LHASET = 29040
458 LHANAME=LHAPATH(1:LHAPATHLEN)//'/MRST98lo.LHgrid'
459 ELSEIF((LHAINPUT .GE. 29050) .AND. (LHAINPUT .LE. 29055)) THEN
460 LHASET = 29050
461 LHANAME=LHAPATH(1:LHAPATHLEN)//'/MRST98nlo.LHgrid'
462 ELSEIF((LHAINPUT .GE. 29060) .AND. (LHAINPUT .LE. 29065)) THEN
463 LHASET = 29060
464 LHANAME=LHAPATH(1:LHAPATHLEN)//'/MRST98dis.LHgrid'
465 ELSEIF((LHAINPUT .GE. 29070) .AND. (LHAINPUT .LE. 29071)) THEN
466 LHASET = 29070
467 LHANAME=LHAPATH(1:LHAPATHLEN)//'/MRST98ht.LHgrid'
468 ELSE
469 WRITE(LHAPRINT,5150) LHASET
470 STOP
471 ENDIF
472C...Fermi Family
473 ELSEIF((LHAINPUT .GE. 30000) .AND. (LHAINPUT .LE. 39999)) THEN
474 IF((LHAINPUT .GE. 30100) .AND. (LHAINPUT .LE. 30200)) THEN
475 LHASET = 30100
476 LHANAME=LHAPATH(1:LHAPATHLEN)//'/Fermi2002_100.LHpdf'
477 ELSEIF((LHAINPUT .GE. 31000) .AND. (LHAINPUT .LE. 32000)) THEN
478 LHASET = 31000
479 LHANAME=LHAPATH(1:LHAPATHLEN)//'/Fermi2002_1000.LHpdf'
480 ELSE
481 WRITE(LHAPRINT,5150) LHASET
482 STOP
483 ENDIF
484C...Alekhin Family
485 ELSEIF((LHAINPUT .GE. 40000) .AND. (LHAINPUT .LE. 49999)) THEN
486 IF((LHAINPUT .GE. 40100) .AND. (LHAINPUT .LE. 40200)) THEN
487 LHASET = 40100
488 LHANAME=LHAPATH(1:LHAPATHLEN)//'/Alekhin_100.LHpdf'
489 ELSEIF((LHAINPUT .GE. 41000) .AND. (LHAINPUT .LE. 41999)) THEN
490 LHASET = 41000
491 LHANAME=LHAPATH(1:LHAPATHLEN)//'/Alekhin_1000.LHpdf'
492 ELSEIF((LHAINPUT .GE. 40350) .AND. (LHAINPUT .LE. 40367)) THEN
493 LHASET = 40350
494 LHANAME=LHAPATH(1:LHAPATHLEN)//'/a02m_lo.LHgrid'
495 XMIN = 1.0D-7
496 Q2MIN = 0.8D0
497 Q2MAX = 2.0D08
498 ELSEIF((LHAINPUT .GE. 40450) .AND. (LHAINPUT .LE. 40467)) THEN
499 LHASET = 40450
500 LHANAME=LHAPATH(1:LHAPATHLEN)//'/a02m_nlo.LHgrid'
501 XMIN = 1.0D-7
502 Q2MIN = 0.8D0
503 Q2MAX = 2.0D08
504 ELSEIF((LHAINPUT .GE. 40550) .AND. (LHAINPUT .LE. 40567)) THEN
505 LHASET = 40550
506 LHANAME=LHAPATH(1:LHAPATHLEN)//'/a02m_nnlo.LHgrid'
507 XMIN = 1.0D-7
508 Q2MIN = 0.8D0
509 Q2MAX = 2.0D08
510 ELSE
511 WRITE(LHAPRINT,5150) LHASET
512 STOP
513 ENDIF
514C...Botje Family
515 ELSEIF((LHAINPUT .GE. 50000) .AND. (LHAINPUT .LE. 59999)) THEN
516 IF((LHAINPUT .GE. 50100) .AND. (LHAINPUT .LE. 50200)) THEN
517 LHASET = 50100
518 LHANAME=LHAPATH(1:LHAPATHLEN)//'/Botje_100.LHpdf'
519 ELSEIF((LHAINPUT .GE. 51000) .AND. (LHAINPUT .LE. 51999)) THEN
520 LHASET = 51000
521 LHANAME=LHAPATH(1:LHAPATHLEN)//'/Botje_1000.LHpdf'
522 ELSE
523 WRITE(LHAPRINT,5150) LHASET
524 STOP
525 ENDIF
526C...ZEUS Family
527 ELSEIF((LHAINPUT .GE. 60000) .AND. (LHAINPUT .LE. 69999)) THEN
528 Q2MIN = 0.3D0
529 Q2MAX = 2.0D05
530 IF((LHAINPUT .GE. 60000) .AND. (LHAINPUT .LE. 60022)) THEN
531 LHASET = 60000
532 LHANAME=LHAPATH(1:LHAPATHLEN)//'/ZEUS2002_TR.LHpdf'
533 ELSEIF((LHAINPUT .GE. 60100) .AND. (LHAINPUT .LE. 60122)) THEN
534 LHASET = 60100
535 LHANAME=LHAPATH(1:LHAPATHLEN)//'/ZEUS2002_ZM.LHpdf'
536 ELSEIF((LHAINPUT .GE. 60200) .AND. (LHAINPUT .LE. 60222)) THEN
537 LHASET = 60200
538 LHANAME=LHAPATH(1:LHAPATHLEN)//'/ZEUS2002_FF.LHpdf'
539 ELSEIF((LHAINPUT .GE. 60300) .AND. (LHAINPUT .LE. 60322)) THEN
540 LHASET = 60300
541 LHANAME=LHAPATH(1:LHAPATHLEN)//'/ZEUS2005_ZJ.LHpdf'
542 ELSE
543 WRITE(LHAPRINT,5150) LHASET
544 STOP
545 ENDIF
546C...H1 Family
547 ELSEIF((LHAINPUT .GE. 70000) .AND. (LHAINPUT .LE. 79999)) THEN
548 Q2MIN = 1.5D0
549 Q2MAX = 3.5D04
550 XMIN = 5.7D-5
551 IF((LHAINPUT .GE. 70050) .AND. (LHAINPUT .LE. 70050)) THEN
552 LHASET = 70050
553 LHANAME=LHAPATH(1:LHAPATHLEN)//'/H12000ms.LHgrid'
554 ELSEIF((LHAINPUT .GE. 70051) .AND. (LHAINPUT .LE. 70070)) THEN
555 LHASET = 70050
556 LHANAME=LHAPATH(1:LHAPATHLEN)//'/H12000msE.LHgrid'
557 ELSEIF((LHAINPUT .GE. 70150) .AND. (LHAINPUT .LE. 70150)) THEN
558 LHASET = 70150
559 LHANAME=LHAPATH(1:LHAPATHLEN)//'/H12000dis.LHgrid'
560 ELSEIF((LHAINPUT .GE. 70151) .AND. (LHAINPUT .LE. 70170)) THEN
561 LHASET = 70150
562 LHANAME=LHAPATH(1:LHAPATHLEN)//'/H12000disE.LHgrid'
563 ELSEIF((LHAINPUT .GE. 70250) .AND. (LHAINPUT .LE. 70250)) THEN
564 LHASET = 70250
565 LHANAME=LHAPATH(1:LHAPATHLEN)//'/H12000lo.LHgrid'
566 ELSEIF((LHAINPUT .GE. 70251) .AND. (LHAINPUT .LE. 70270)) THEN
567 LHASET = 70250
568 LHANAME=LHAPATH(1:LHAPATHLEN)//'/H12000loE.LHgrid'
569c tempoararily removed on returning to original H!2000 files
570c ELSEIF((LHAINPUT .GE. 70350) .AND. (LHAINPUT .LE. 70350)) THEN
571c LHASET = 70350
572c LHANAME=LHAPATH(1:LHAPATHLEN)//'/H12000lo2.LHgrid'
573c ELSEIF((LHAINPUT .GE. 70351) .AND. (LHAINPUT .LE. 70370)) THEN
574c LHASET = 70350
575c LHANAME=LHAPATH(1:LHAPATHLEN)//'/H12000lo2E.LHgrid'
576 ELSE
577 WRITE(LHAPRINT,5150) LHASET
578 STOP
579 ENDIF
580C...GRV Family
581 ELSEIF((LHAINPUT .GE. 80000) .AND. (LHAINPUT .LE. 89999)) THEN
582 Q2MIN = 0.8D0
583 Q2MAX = 2.0D06
584 XMIN = 1.0D-9
585 IF((LHAINPUT .GE. 80050) .AND. (LHAINPUT .LE. 80051)) THEN
586 LHASET = 80050
587 LHANAME=LHAPATH(1:LHAPATHLEN)//'/GRV98nlo.LHgrid'
588 ELSEIF((LHAINPUT .GE. 80060) .AND. (LHAINPUT .LE. 80060)) THEN
589 LHASET = 80060
590 LHANAME=LHAPATH(1:LHAPATHLEN)//'/GRV98lo.LHgrid'
591 ELSE
592 WRITE(LHAPRINT,5150) LHASET
593 STOP
594 ENDIF
595C...
596C...Pions
597C...
598C...OW-PI Family
599 ELSEIF((LHAINPUT .GE. 210) .AND. (LHAINPUT .LE. 212)) THEN
600 Q2MIN = 4.0D0
601 Q2MAX = 2.0D03
602 XMIN = 5.0D-03
603 XMAX = 0.9998D0
604 IF((LHAINPUT .GE. 210) .AND. (LHAINPUT .LE. 212)) THEN
605 LHASET = 210
606 LHANAME=LHAPATH(1:LHAPATHLEN)//'/OWPI.LHgrid'
607 ELSE
608 WRITE(LHAPRINT,5150) LHASET
609 STOP
610 ENDIF
611C...SMRS-PI Family
612 ELSEIF((LHAINPUT .GE. 230) .AND. (LHAINPUT .LE. 233)) THEN
613 Q2MIN = 5.0D0
614 Q2MAX = 1.31D06
615 XMIN = 1.0D-05
616 XMAX = 0.9998D0
617 IF((LHAINPUT .GE. 230) .AND. (LHAINPUT .LE. 233)) THEN
618 LHASET = 230
619 LHANAME=LHAPATH(1:LHAPATHLEN)//'/SMRSPI.LHgrid'
620 ELSE
621 WRITE(LHAPRINT,5150) LHASET
622 STOP
623 ENDIF
624C...GRV-PI Family
625 ELSEIF((LHAINPUT .GE. 250) .AND. (LHAINPUT .LE. 252)) THEN
626 Q2MAX = 1.00D06
627 XMIN = 1.0D-05
628 XMAX = 0.9998D0
629 IF((LHAINPUT .GE. 250) .AND. (LHAINPUT .LE. 251)) THEN
630 Q2MIN = 3.0D-1
631 LHASET = 250
632 LHANAME=LHAPATH(1:LHAPATHLEN)//'/GRVPI1.LHgrid'
633 ELSEIF((LHAINPUT .GE. 252) .AND. (LHAINPUT .LE. 252)) THEN
634 Q2MIN = 2.5D-1
635 LHASET = 252
636 LHANAME=LHAPATH(1:LHAPATHLEN)//'/GRVPI0.LHgrid'
637 ELSE
638 WRITE(LHAPRINT,5150) LHASET
639 STOP
640 ENDIF
641C...ABFKW-PI Family
642 ELSEIF((LHAINPUT .GE. 260) .AND. (LHAINPUT .LE. 263)) THEN
643 Q2MIN = 2.0D0
644 Q2MAX = 1.00D08
645 XMIN = 1.0D-03
646 XMAX = 0.9998D0
647 IF((LHAINPUT .GE. 260) .AND. (LHAINPUT .LE. 263)) THEN
648 LHASET = 260
649 LHANAME=LHAPATH(1:LHAPATHLEN)//'/ABFKWPI.LHgrid'
650 ELSE
651 WRITE(LHAPRINT,5150) LHASET
652 STOP
653 ENDIF
654C...
655C...Photons
656C...
657C...DO-G Family
658 ELSEIF((LHAINPUT .GE. 310) .AND. (LHAINPUT .LE. 312)) THEN
659 Q2MIN = 1.0D01
660 Q2MAX = 1.00D04
661 XMIN = 1.0D-05
662 XMAX = 0.9D0
663 IF((LHAINPUT .GE. 310) .AND. (LHAINPUT .LE. 311)) THEN
664 LHASET = 310
665 LHANAME=LHAPATH(1:LHAPATHLEN)//'/DOG0.LHgrid'
666 ELSEIF((LHAINPUT .GE. 312) .AND. (LHAINPUT .LE. 312)) THEN
667 LHASET = 312
668 LHANAME=LHAPATH(1:LHAPATHLEN)//'/DOG1.LHgrid'
669 ELSE
670 WRITE(LHAPRINT,5150) LHASET
671 STOP
672 ENDIF
673C...DG-G Family
674 ELSEIF((LHAINPUT .GE. 320) .AND. (LHAINPUT .LE. 324)) THEN
675 XMIN = 1.0D-05
676 XMAX = 0.9998D0
677 LHASET = 320
678 IF((LHAINPUT .GE. 320) .AND. (LHAINPUT .LE. 321)) THEN
679 Q2MIN = 1.0D0
680 Q2MAX = 1.0D04
681c LHASET = 320
682 LHANAME=LHAPATH(1:LHAPATHLEN)//'/DGG.LHgrid'
683 ELSEIF((LHAINPUT .GE. 322) .AND. (LHAINPUT .LE. 322)) THEN
684 Q2MIN = 1.0D0
685 Q2MAX = 5.0D01
686c LHASET = 322
687 LHANAME=LHAPATH(1:LHAPATHLEN)//'/DGG.LHgrid'
688 ELSEIF((LHAINPUT .GE. 323) .AND. (LHAINPUT .LE. 323)) THEN
689 Q2MIN = 2.0D1
690 Q2MAX = 5.0D02
691c LHASET = 323
692 LHANAME=LHAPATH(1:LHAPATHLEN)//'/DGG.LHgrid'
693 ELSEIF((LHAINPUT .GE. 324) .AND. (LHAINPUT .LE. 324)) THEN
694 Q2MIN = 2.0D2
695 Q2MAX = 1.0D04
696c LHASET = 324
697 LHANAME=LHAPATH(1:LHAPATHLEN)//'/DGG.LHgrid'
698 ELSE
699 WRITE(LHAPRINT,5150) LHASET
700 STOP
701 ENDIF
702C...LAC/GAL-G Family
703 ELSEIF((LHAINPUT .GE. 330) .AND. (LHAINPUT .LE. 334)) THEN
704 Q2MIN = 4.0D00
705 Q2MAX = 1.0D05
706 XMIN = 1.0D-04
707 XMAX = 0.9998D0
708 LHASET = 330
709 IF((LHAINPUT .GE. 330) .AND. (LHAINPUT .LE. 332)) THEN
710c LHASET = 330
711 LHANAME=LHAPATH(1:LHAPATHLEN)//'/LACG.LHgrid'
712 ELSEIF((LHAINPUT .GE. 333) .AND. (LHAINPUT .LE. 333)) THEN
713 Q2MIN = 1.0D00
714c LHASET = 333
715 LHANAME=LHAPATH(1:LHAPATHLEN)//'/LACG.LHgrid'
716 ELSEIF((LHAINPUT .GE. 334) .AND. (LHAINPUT .LE. 334)) THEN
717 Q2MIN = 4.0D00
718c LHASET = 334
719 LHANAME=LHAPATH(1:LHAPATHLEN)//'/LACG.LHgrid'
720 ELSE
721 WRITE(LHAPRINT,5150) LHASET
722 STOP
723 ENDIF
724C...GSG/GSG96-G Family
725 ELSEIF((LHAINPUT .GE. 340) .AND. (LHAINPUT .LE. 345)) THEN
726 Q2MIN = 5.3D00
727 Q2MAX = 1.0D08
728 XMIN = 5.0D-04
729 XMAX = 0.9998D0
730 IF((LHAINPUT .GE. 340) .AND. (LHAINPUT .LE. 341)) THEN
731 LHASET = 340
732 LHANAME=LHAPATH(1:LHAPATHLEN)//'/GSG1.LHgrid'
733 ELSEIF((LHAINPUT .GE. 342) .AND. (LHAINPUT .LE. 343)) THEN
734 LHASET = 342
735 LHANAME=LHAPATH(1:LHAPATHLEN)//'/GSG0.LHgrid'
736 ELSEIF((LHAINPUT .GE. 344) .AND. (LHAINPUT .LE. 344)) THEN
737 LHASET = 344
738 LHANAME=LHAPATH(1:LHAPATHLEN)//'/GSG961.LHgrid'
739 ELSEIF((LHAINPUT .GE. 345) .AND. (LHAINPUT .LE. 345)) THEN
740 LHASET = 345
741 LHANAME=LHAPATH(1:LHAPATHLEN)//'/GSG960.LHgrid'
742 ELSE
743 WRITE(LHAPRINT,5150) LHASET
744 STOP
745 ENDIF
746C...GRV-G Family
747 ELSEIF((LHAINPUT .GE. 350) .AND. (LHAINPUT .LE. 354)) THEN
748 Q2MIN = 3.0D-1
749 Q2MAX = 1.0D06
750 XMIN = 1.0D-05
751 XMAX = 0.9998D0
752 IF((LHAINPUT .GE. 350) .AND. (LHAINPUT .LE. 352)) THEN
753 LHASET = 350
754 LHANAME=LHAPATH(1:LHAPATHLEN)//'/GRVG1.LHgrid'
755 ELSEIF((LHAINPUT .GE. 353) .AND. (LHAINPUT .LE. 353)) THEN
756 Q2MIN = 2.5D-1
757 LHASET = 353
758 LHANAME=LHAPATH(1:LHAPATHLEN)//'/GRVG0.LHgrid'
759 ELSEIF((LHAINPUT .GE. 354) .AND. (LHAINPUT .LE. 354)) THEN
760 Q2MIN = 6.0D-1
761 Q2MAX = 5.0D04
762 LHASET = 354
763 LHANAME=LHAPATH(1:LHAPATHLEN)//'/GRVG0.LHgrid'
764 ELSE
765 WRITE(LHAPRINT,5150) LHASET
766 STOP
767 ENDIF
768C...ACFGP-G Family
769 ELSEIF((LHAINPUT .GE. 360) .AND. (LHAINPUT .LE. 363)) THEN
770 Q2MIN = 2.0D00
771 Q2MAX = 5.5D05
772 XMIN = 1.37D-03
773 XMAX = 0.9998D0
774 IF((LHAINPUT .GE. 360) .AND. (LHAINPUT .LE. 363)) THEN
775 LHASET = 360
776 LHANAME=LHAPATH(1:LHAPATHLEN)//'/ACFGPG.LHgrid'
777 ELSE
778 WRITE(LHAPRINT,5150) LHASET
779 STOP
780 ENDIF
781C...WHIT-G Family
782 ELSEIF((LHAINPUT .GE. 380) .AND. (LHAINPUT .LE. 386)) THEN
783 Q2MIN = 4.0D00
784 Q2MAX = 2.5D03
785 XMIN = 1.0D-03
786 XMAX = 0.9998D0
787 IF((LHAINPUT .GE. 380) .AND. (LHAINPUT .LE. 386)) THEN
788 LHASET = 380
789 LHANAME=LHAPATH(1:LHAPATHLEN)//'/WHITG.LHgrid'
790 ELSE
791 WRITE(LHAPRINT,5150) LHASET
792 STOP
793 ENDIF
794C...SAS-G Family
795 ELSEIF((LHAINPUT .GE. 390) .AND. (LHAINPUT .LE. 398)) THEN
796 Q2MAX = 5.0D04
797 XMIN = 1.0D-05
798 XMAX = 0.9998D0
799 LHASET = 390
800 IF((LHAINPUT .GE. 390) .AND. (LHAINPUT .LE. 392)) THEN
801 Q2MIN = 3.6D-1
802c LHASET = 390
803 LHANAME=LHAPATH(1:LHAPATHLEN)//'/SASG.LHgrid'
804 ELSEIF((LHAINPUT .GE. 393) .AND. (LHAINPUT .LE. 394)) THEN
805 Q2MIN = 4.0D00
806c LHASET = 393
807 LHANAME=LHAPATH(1:LHAPATHLEN)//'/SASG.LHgrid'
808 ELSEIF((LHAINPUT .GE. 395) .AND. (LHAINPUT .LE. 396)) THEN
809 Q2MIN = 3.6D-1
810c LHASET = 395
811 LHANAME=LHAPATH(1:LHAPATHLEN)//'/SASG.LHgrid'
812 ELSEIF((LHAINPUT .GE. 397) .AND. (LHAINPUT .LE. 398)) THEN
813 Q2MIN = 4.0D00
814c LHASET = 397
815 LHANAME=LHAPATH(1:LHAPATHLEN)//'/SASG.LHgrid'
816 ELSE
817 WRITE(LHAPRINT,5150) LHASET
818 STOP
819 ENDIF
820C...Unknown Family ?! Giving up
821 ELSE
822 WRITE(LHAPRINT,5150) LHASET
823 STOP
824 ENDIF
825
826 LHAMEMB=LHAINPUT-LHASET
827c....Now work out if we have already called this set/member
828 iset = 0
829 do j=1,nsets
830 if(lhaname.eq.lhanames(j).and.
831 + lhamemb.eq.lhamembers(j)) then
832 iset = j
833 endif
834 enddo
835 if(iset.eq.0) then
836 nsets=nsets+1
837 if(nsets.gt.nmxset) then
838 if(LHASILENT.ne.1) then
839 print *,'WARNING:too many sets initialised'
840 print *,'overwriting from set 1 again'
841 endif
842 nsets = 1
843c stop
844 endif
845 iset=nsets
846 lhanames(iset)=lhaname
847 lhanumbers(iset)=lhainput
848 lhamembers(iset)=lhamemb
849 xxmin(iset)=xmin
850 xxmax(iset)=xmax
851 qq2min(iset)=q2min
852 qq2max(iset)=q2max
853 CALL INITPDFSETM(iset,LHANAME)
854 CALL NUMBERPDFM(iset,LHAALLMEM)
855 IF(LHASILENT .NE. 1) THEN
856 WRITE(LHAPRINT,5151)
857 WRITE(LHAPRINT,5152) LHANAME
858 WRITE(LHAPRINT,5153) LHAALLMEM
859 WRITE(LHAPRINT,5154)
860 ENDIF
861 IF ((LHAMEMB.LT.0) .OR. (LHAMEMB.GT.LHAALLMEM)) THEN
862 WRITE(LHAPRINT,5155) LHAMEMB
863 WRITE(LHAPRINT,5156) LHAALLMEM
864 STOP
865 ENDIF
866
867c print *,'calling initpdf',lhamemb
868c print *,'calling initpdfm ',iset,lhaname,lhamemb
869c print *,'LHAGLUE .... initializing set,member ',iset,lhamemb
870 CALL INITPDFM(iset,LHAMEMB)
871 endif
872c... the rest is done every time pdfset is called
873c print *,'setting nset to:',iset
874 call setnset(iset)
875 call setnmem(iset,lhamemb)
876 xmin = xxmin(iset)
877 xmax = xxmax(iset)
878 q2min=qq2min(iset)
879 q2max=qq2max(iset)
880 call GetLam4M(iset,LHAMEMB,qcdl4)
881 call GetLam5M(iset,LHAMEMB,qcdl5)
882
883 QMZ = 91.1876D0
884 alphasLHA = alphasPDFM(iset,QMZ)
885 IF(LHASILENT .NE. 1)
886 > WRITE(LHAPRINT,5158) alphasLHA
887
888 IF(LHAPARM(17).EQ.'LHAPDF') THEN
889 NPTYPEPDFL = 1 ! Proton PDFs
890 NFLPDFL = 4
891 QCDLHA4 = QCDL4
892 QCDLHA5 = QCDL5
893 IF(LHASILENT .NE. 1)
894 > WRITE(LHAPRINT,5159) QCDL4, QCDL5
895 ELSE
896 NPTYPEPDFL = 1 ! Proton PDFs
897 NFLPDFL = 4
898 ALAMBDA = 0.192D0
899 QCDLHA4 = ALAMBDA
900 QCDLHA5 = ALAMBDA
901 IF(PARM(1).EQ.'NPTYPE') THEN ! PYTHIA
902 QCDL4 = ALAMBDA
903 QCDL5 = ALAMBDA
904 ENDIF
905 ENDIF
906C...Formats for initialization information.
907 5150 FORMAT(1X,'WRONG LHAPDF set number =',I12,' given! STOP EXE!')
908 5151 FORMAT(1X,'==============================================')
909 5152 FORMAT(1X,'PDFset name ',A80)
910 5153 FORMAT(1X,'with ',I10,' members')
911 5154 FORMAT(1X,'==== initialized. ===========================')
912 5155 FORMAT(1X,'LHAPDF problem => YOU asked for member = ',I10)
913 5156 FORMAT(1X,'Valid range is: 0 - ',I10,' Execution stopped.')
914 5157 FORMAT(1X,'Number of flavors for PDF is:',I4)
915 5158 FORMAT(1X,'Strong coupling at Mz for PDF is:',F9.5)
916 5159 FORMAT(1X,'Will use for PYTHIA QCDL4, QCDL5:',2F9.5)
917
918 RETURN
919 END
920
921c********************************************************************
922c -- STRUCTA
923c -- copy of PDFLIB to use the eks98 nuclear correction factors
924
925 SUBROUTINE STRUCTA(X,Q,A,UPV,DNV,USEA,DSEA,STR,CHM,BOT,TOP,GLU)
926 IMPLICIT DOUBLE PRECISION (A-H,O-Z)
927 CALL EKS98(X,Q,A,RUV,RDV,RU,RD,RS,RC,RB,RT,RG)
928 CALL STRUCTM(X,Q,UPV,DNV,USEA,DSEA,STR,CHM,BOT,TOP,GLU)
929 UPV = RUV*UPV
930 DNV = RDV*DNV
931 USEA = RU*USEA
932 DSEA = RD*DSEA
933 STR = RS*STR
934 CHM = RC*CHM
935 BOT = RB*BOT
936 TOP = RT*TOP
937 GLU = RG*GLU
938 RETURN
939 END
940
941C*********************************************************************
942
943C...STRUCTM
944C...Gives parton distributions according to the LHAPDF interface.
945C...Two evolution codes used:
946C... EVLCTEQ for CTEQ PDF sets
947C... QCDNUM for Other PDF sets
948C...
949C...Author: Dimitri Bourilkov bourilkov@mailaps.org
950C...
951C...v4.0 21-Mar-2005 Photon/pion/new p PDFs, updated for LHAPDF v4
952C...v3.0 23-Jan-2004
953C...
954C...interface to LHAPDF library
955
956 SUBROUTINE STRUCTM(X,Q,UPV,DNV,USEA,DSEA,STR,CHM,BOT,TOP,GLU)
957
958C...Double precision and integer declarations.
959 IMPLICIT DOUBLE PRECISION(A-H, O-Z)
960 IMPLICIT INTEGER(I-N)
961C include 'parmsetup.inc'
962 integer noemax,nopmax,npfmax,nofmax,linemax,nmxset
963 parameter (noemax=1000,nopmax=40,npfmax=10,nofmax=10,linemax=20)
964c nmxset is the max number of sets that can be initialised at one time ---- added V5
965 parameter (nmxset=3)
966C...Commonblocks.
967 COMMON/PYDAT1/MSTU(200),PARU(200),MSTJ(200),PARJ(200)
968 SAVE /PYDAT1/
969C...Interface to LHAPDFLIB.
970C include 'pathsetup.inc'
971 character*232 lhapath
972 character*272 lhaname,lhanames(nmxset)
973c CHARACTER*172 LHANAME
974 INTEGER LHASET, LHAMEMB
975 COMMON/LHAPDF/LHANAME, LHASET, LHAMEMB
976 SAVE /LHAPDF/
977c...added next 2 lines for structp fix
978 integer LHAMEMBERS(nmxset),LHANUMBERS(nmxset)
979 common/LHASETS/LHANAMES,LHANUMBERS,LHAMEMBERS,nsets
980c
981 DOUBLE PRECISION QCDLHA4, QCDLHA5
982 INTEGER NFLLHA
983 COMMON/LHAPDFR/QCDLHA4, QCDLHA5, NFLLHA
984 SAVE /LHAPDFR/
985 CHARACTER*20 LHAPARM(20)
986 DOUBLE PRECISION LHAVALUE(20)
987 COMMON/LHACONTROL/LHAPARM,LHAVALUE
988 SAVE/LHACONTROL/
989 INTEGER LHAEXTRP
990 COMMON/LHAPDFE/LHAEXTRP
991 SAVE /LHAPDFE/
992 DOUBLE PRECISION XMINNUM,XMAXNUM,Q2MINNUM,Q2MAXNUM,TOTNUM,
993 > XMINNUP,XMAXNUP,Q2MINNUP,Q2MAXNUP,TOTNUP
994 COMMON/LHAGLSTA/ XMINNUM,XMAXNUM,Q2MINNUM,Q2MAXNUM,TOTNUM,
995 > XMINNUP,XMAXNUP,Q2MINNUP,Q2MAXNUP,TOTNUP
996 SAVE/LHAGLSTA/
997C...Interface to PDFLIB.
998 COMMON/W50513/XMIN,XMAX,Q2MIN,Q2MAX
999 SAVE /W50513/
1000 DOUBLE PRECISION XMIN,XMAX,Q2MIN,Q2MAX
1001C...Local variables
1002 DOUBLE PRECISION UPV,DNV,USEA,DSEA,STR,CHM,BOT,TOP,GLU
1003 DOUBLE PRECISION X,Q,F(-6:6)
1004
1005 Q2 = Q**2
1006C...Statistics
1007 IF(LHAPARM(16).NE.'NOSTAT') THEN
1008 TOTNUM = TOTNUM+1.D0
1009 IF(X .LT. XMIN) XMINNUM = XMINNUM+1.D0
1010 IF(X .GT. XMAX) XMAXNUM = XMAXNUM+1.D0
1011 IF(Q2 .LT. Q2MIN) Q2MINNUM = Q2MINNUM+1.D0
1012 IF(Q2 .GT. Q2MAX) Q2MAXNUM = Q2MAXNUM+1.D0
1013 ENDIF
1014
1015C...Range of validity e.g. 10^-6 < x < 1, Q2MIN < Q^2 extended by
1016C...freezing x*f(x,Q2) at borders.
1017 IF(LHAEXTRP .NE. 1) THEN ! safe mode == "freeze"
1018 XIN=MAX(XMIN,MIN(XMAX,X))
1019 Q=SQRT(MAX(0D0,Q2MIN,MIN(Q2MAX,Q2)))
1020 ELSE ! adventurous mode == OWN RISK !
1021 XIN=X
1022 ENDIF
1023
1024 call getnset(iset)
1025c print *,'calling evolvepdfm:',iset
1026C
1027C...fix to allow STRUCTM to work for photon PDFs (Herwig does this)
1028C...set P2 = 0.0d0 and IP2 = 0
1029 if(LHANUMBERS(iset).ge.300.and.LHANUMBERS(iset).le.399) then
1030 P2 = 0.0d0
1031 IP2 = 0
1032 CALL EVOLVEPDFPM(iset,XIN,Q,P2,IP2,F)
1033 else
1034 CALL EVOLVEPDFM(iset,XIN,Q,F)
1035 endif
1036 GLU = F(0)
1037 DSEA = F(-1)
1038 DNV = F(1) - DSEA
1039 USEA = F(-2)
1040 UPV = F(2) - USEA
1041 STR = F(3)
1042 CHM = F(4)
1043 BOT = F(5)
1044 TOP = F(6)
1045
1046 RETURN
1047 END
1048
1049C*********************************************************************
1050
1051C...STRUCTP
1052C...Gives parton distributions according to the LHAPDF interface.
1053C...Used for photons.
1054C...
1055C...v4.0 21-Mar-2005 Photon/pion/new p PDFs, updated for LHAPDF v4
1056C...
1057C...interface to LHAPDF library
1058
1059 SUBROUTINE STRUCTP
1060 > (X,Q2,P2,IP2,UPV,DNV,USEA,DSEA,STR,CHM,BOT,TOP,GLU)
1061
1062C...Double precision and integer declarations.
1063 IMPLICIT DOUBLE PRECISION(A-H, O-Z)
1064 IMPLICIT INTEGER(I-N)
1065C include 'parmsetup.inc'
1066 integer noemax,nopmax,npfmax,nofmax,linemax,nmxset
1067 parameter (noemax=1000,nopmax=40,npfmax=10,nofmax=10,linemax=20)
1068c nmxset is the max number of sets that can be initialised at one time ---- added V5
1069 parameter (nmxset=3)
1070C...Commonblocks.
1071 COMMON/PYDAT1/MSTU(200),PARU(200),MSTJ(200),PARJ(200)
1072 SAVE /PYDAT1/
1073C...Interface to LHAPDFLIB.
1074C include 'pathsetup.inc'
1075 character*232 lhapath
1076 character*272 lhaname,lhanames(nmxset)
1077c CHARACTER*172 LHANAME
1078 INTEGER LHASET, LHAMEMB
1079 COMMON/LHAPDF/LHANAME, LHASET, LHAMEMB
1080 SAVE /LHAPDF/
1081 DOUBLE PRECISION QCDLHA4, QCDLHA5
1082 INTEGER NFLLHA
1083 COMMON/LHAPDFR/QCDLHA4, QCDLHA5, NFLLHA
1084 SAVE /LHAPDFR/
1085 CHARACTER*20 LHAPARM(20)
1086 DOUBLE PRECISION LHAVALUE(20)
1087 COMMON/LHACONTROL/LHAPARM,LHAVALUE
1088 SAVE/LHACONTROL/
1089 INTEGER LHAEXTRP
1090 COMMON/LHAPDFE/LHAEXTRP
1091 SAVE /LHAPDFE/
1092 DOUBLE PRECISION XMINNUM,XMAXNUM,Q2MINNUM,Q2MAXNUM,TOTNUM,
1093 > XMINNUP,XMAXNUP,Q2MINNUP,Q2MAXNUP,TOTNUP
1094 COMMON/LHAGLSTA/ XMINNUM,XMAXNUM,Q2MINNUM,Q2MAXNUM,TOTNUM,
1095 > XMINNUP,XMAXNUP,Q2MINNUP,Q2MAXNUP,TOTNUP
1096 SAVE/LHAGLSTA/
1097C...Interface to PDFLIB.
1098 COMMON/W50513/XMIN,XMAX,Q2MIN,Q2MAX
1099 SAVE /W50513/
1100 DOUBLE PRECISION XMIN,XMAX,Q2MIN,Q2MAX
1101C...Local variables
1102 DOUBLE PRECISION UPV,DNV,USEA,DSEA,STR,CHM,BOT,TOP,GLU
1103 DOUBLE PRECISION X,Q,F(-6:6)
1104
1105C...Statistics
1106 IF(LHAPARM(16).NE.'NOSTAT') THEN
1107 TOTNUP = TOTNUP+1.D0
1108 IF(X .LT. XMIN) XMINNUP = XMINNUP+1.D0
1109 IF(X .GT. XMAX) XMAXNUP = XMAXNUP+1.D0
1110 IF(Q2 .LT. Q2MIN) Q2MINNUP = Q2MINNUP+1.D0
1111 IF(Q2 .GT. Q2MAX) Q2MAXNUP = Q2MAXNUP+1.D0
1112 ENDIF
1113
1114C...Range of validity e.g. 10^-6 < x < 1, Q2MIN < Q^2 extended by
1115C...freezing x*f(x,Q2) at borders.
1116 Q = DSQRT(Q2)
1117 IF(LHAEXTRP .NE. 1) THEN ! safe mode == "freeze"
1118 XIN=MAX(XMIN,MIN(XMAX,X))
1119 Q=SQRT(MAX(0D0,Q2MIN,MIN(Q2MAX,Q2)))
1120 ELSE ! adventurous mode == OWN RISK !
1121 XIN=X
1122 ENDIF
1123 call getnset(iset)
1124 CALL EVOLVEPDFPM(iset,XIN,Q,P2,IP2,F)
1125
1126 GLU = F(0)
1127 DSEA = F(-1)
1128 DNV = F(1) - DSEA
1129 USEA = F(-2)
1130 UPV = F(2) - USEA
1131 STR = F(3)
1132 CHM = F(4)
1133 BOT = F(5)
1134 TOP = F(6)
1135
1136 RETURN
1137 END
1138
1139C*********************************************************************
1140
1141C...PDFSTA
1142C...For statistics ON structure functions (under/over-flows)
1143C...
1144C...Author: Dimitri Bourilkov bourilkov@mailaps.org
1145C...
1146C...
1147C...first introduced in v4.0 28-Apr-2005
1148C...
1149
1150 SUBROUTINE PDFSTA
1151
1152C...Double precision and integer declarations.
1153 IMPLICIT DOUBLE PRECISION(A-H, O-Z)
1154 IMPLICIT INTEGER(I-N)
1155C...Interface to LHAPDFLIB.
1156 DOUBLE PRECISION XMINNUM,XMAXNUM,Q2MINNUM,Q2MAXNUM,TOTNUM,
1157 > XMINNUP,XMAXNUP,Q2MINNUP,Q2MAXNUP,TOTNUP
1158 COMMON/LHAGLSTA/ XMINNUM,XMAXNUM,Q2MINNUM,Q2MAXNUM,TOTNUM,
1159 > XMINNUP,XMAXNUP,Q2MINNUP,Q2MAXNUP,TOTNUP
1160 SAVE/LHAGLSTA/
1161
1162 PRINT *
1163 PRINT *,'===== PDFSTA statistics for PDF under/over-flows ===='
1164 PRINT *
1165 PRINT *,'====== STRUCTM statistics for nucleon/pion PDFs ====='
1166 PRINT *
1167 PRINT *,' total # of calls ',TOTNUM
1168 IF(TOTNUM .GT. 0.D0) THEN
1169 PERCBELOW = 100.D0*XMINNUM/TOTNUM
1170 PERCABOVE = 100.D0*XMAXNUM/TOTNUM
1171 PRINT *,' X below PDF min ',XMINNUM,' or ',PERCBELOW, ' %'
1172 PRINT *,' X above PDF max ',XMAXNUM,' or ',PERCABOVE, ' %'
1173 PERCBELOW = 100.D0*Q2MINNUM/TOTNUM
1174 PERCABOVE = 100.D0*Q2MAXNUM/TOTNUM
1175 PRINT *,' Q2 below PDF min ',Q2MINNUM,' or ',PERCBELOW, ' %'
1176 PRINT *,' Q2 above PDF max ',Q2MAXNUM,' or ',PERCABOVE, ' %'
1177 ENDIF
1178 PRINT *
1179 PRINT *,'========= STRUCTP statistics for photon PDFs ========'
1180 PRINT *
1181 PRINT *,' total # of calls ',TOTNUP
1182 IF(TOTNUP .GT. 0.D0) THEN
1183 PERCBELOW = 100.D0*XMINNUP/TOTNUP
1184 PERCABOVE = 100.D0*XMAXNUP/TOTNUP
1185 PRINT *,' X below PDF min ',XMINNUP,' or ',PERCBELOW, ' %'
1186 PRINT *,' X above PDF max ',XMAXNUP,' or ',PERCABOVE, ' %'
1187 PERCBELOW = 100.D0*Q2MINNUP/TOTNUP
1188 PERCABOVE = 100.D0*Q2MAXNUP/TOTNUP
1189 PRINT *,' Q2 below PDF min ',Q2MINNUP,' or ',PERCBELOW, ' %'
1190 PRINT *,' Q2 above PDF max ',Q2MAXNUP,' or ',PERCABOVE, ' %'
1191 ENDIF
1192 PRINT *
1193
1194 RETURN
1195 END
1196**********************************************************************
1197*
1198* $Id$
1199*
1200* $Log$
1201* Revision 1.3 2006/11/02 10:22:57 hristov
1202* Extracting the BLOCK DATA in a separate file. Changes to make it working on macosx
1203*
1204* Revision 1.2 2006/11/01 12:25:47 hristov
1205* Using LHAPDF instead of PDF
1206*
1207* Revision 1.1 2006/08/07 09:09:40 morsch
1208* LHAPDF 5.2.2 source code.
1209*
1210* Revision 1.7 2005/12/02 14:50:54 whalley
1211* Changes for new CTEQ code/AB sets
1212*
1213* Revision 1.6 2005/10/18 15:35:48 whalley
1214* fix to allow LHAPATH to be user defined as well as lhapdf-config
1215*
1216* Revision 1.5 2005/10/18 11:47:48 whalley
1217* Change to only set LHAPATH once per run
1218*
1219* Revision 1.1.1.2 1996/10/30 08:29:06 cernlib
1220* Version 7.04
1221*
1222* Revision 1.1.1.1 1996/04/12 15:29:26 plothow
1223* Version 7.01
1224*
1225*
1226 SUBROUTINE PFTOPDG(DX,DSCALE,DXPDF)
1227C
1228Cinclude "pdf/expdp.inc"
1229 double precision
1230 + DX,DSCALE,DUPV,DDNV,DUSEA,DDSEA,DSTR,DCHM,DBOT,DTOP,DGL,
1231 + DXPDF(-6:6)
1232C... call STRUCTM in PDFLIB to get flavour content
1233 CALL STRUCTM(DX,DSCALE,
1234 + DUPV,DDNV,DUSEA,DDSEA,DSTR,DCHM,DBOT,DTOP,DGL)
1235C... convert flavour convention of PDFLIB to PDG convention
1236 DXPDF(0) = DGL
1237 DXPDF(1) = DDNV + DDSEA
1238 DXPDF(2) = DUPV + DUSEA
1239 DXPDF(3) = DSTR
1240 DXPDF(4) = DCHM
1241 DXPDF(5) = DBOT
1242 DXPDF(6) = DTOP
1243 DXPDF(-1) = DDSEA
1244 DXPDF(-2) = DUSEA
1245 DXPDF(-3) = DSTR
1246 DXPDF(-4) = DCHM
1247 DXPDF(-5) = DBOT
1248 DXPDF(-6) = DTOP
1249C
1250 RETURN
1251 END
1252****************************************************************************
1253 subroutine setPDFpath(pathname)
1254 implicit real*8 (A-H,O-Z)
1255C include 'parmsetup.inc'
1256 integer noemax,nopmax,npfmax,nofmax,linemax,nmxset
1257 parameter (noemax=1000,nopmax=40,npfmax=10,nofmax=10,linemax=20)
1258c nmxset is the max number of sets that can be initialised at one time ---- added V5
1259 parameter (nmxset=3)
1260 character*(*) pathname
1261C include 'pathsetup.inc'
1262 character*232 lhapath
1263 character*272 lhaname,lhanames(nmxset)
1264c character*132 lhapath
1265 common/LHAPDFC/lhapath
1266 character*20 lhaparm(20)
1267 real*8 lhavalue(20)
1268 common/LHACONTROL/lhaparm,lhavalue
1269 lhaparm(20) = 'LHAPATH'
1270 do j=1,lnblnk(lhapath)
1271 lhapath(j:j)=''
1272 enddo
1273 lhapath = pathname
1274 return
1275 end
1276***********************************************************************
1277 subroutine lhaset(lhaparm2,lhavalue2)
1278 implicit real*8 (a-h,o-z)
1279 character*20 lhaparm(20),lhaparm2(20)
1280 real*8 lhavalue(20),lhavalue2(20)
1281 common/LHACONTROL/lhaparm,lhavalue
1282 do j=1,20
1283 lhaparm(j)=lhaparm2(j)
1284 lhavalue(j)=lhavalue2(j)
1285 enddo
1286 return
1287 end
1288******************************************************************
1289 subroutine setlhaparm(lparm)
1290 implicit real*8 (a-h,o-z)
1291 character*(*) lparm
1292 character*20 lhaparm(20)
1293 real*8 lhavalue(20)
1294 common/LHACONTROL/lhaparm,lhavalue
1295 if(lparm.eq.'NOSTAT') then
1296 lhaparm(16)='NOSTAT'
1297 else if (lparm.eq.'16') then
1298 lhaparm(16)=''
1299 else if (lparm.eq.'LHAPDF') then
1300 lhaparm(17)='LHAPDF'
1301 else if (lparm.eq.'17') then
1302 lhaparm(17)=''
1303 else if (lparm.eq.'EXTRAPOLATE') then
1304 lhaparm(18)='EXTRAPOLATE'
1305 else if (lparm.eq.'18') then
1306 lhaparm(18)=''
1307 else if (lparm.eq.'SILENT') then
1308 lhaparm(19)='SILENT'
1309 else if (lparm.eq.'LOWKEY') then
1310 lhaparm(19)='LOWKEY'
1311 else if (lparm.eq.'19') then
1312 lhaparm(19)=''
1313 else
1314 print *,'WARNING from SetLHAPARM - value',lparm,'
1315 & not recognized!'
1316 endif
1317 return
1318 end
1319***************************************************************