EIC Software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
grv98pa.F
Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file grv98pa.F
1 
2 
3 * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
4 * *
5 * G R V - P R O T O N - P A R A M E T R I Z A T I O N S *
6 * *
7 * 1998 UPDATE *
8 * *
9 * For a detailed explanation see *
10 * M. Glueck, E. Reya, A. Vogt : *
11 * hep-ph/9806404 = DO-TH 98/07 = WUE-ITP-98-019 *
12 * (To appear in Eur. Phys. J. C) *
13 * *
14 * This package contains subroutines returning the light-parton *
15 * distributions in NLO (for the MSbar and DIS schemes) and LO; *
16 * the respective light-parton, charm, and bottom contributions *
17 * to F2(electromagnetic); and the scale dependence of alpha_s. *
18 * *
19 * The parton densities and F2 values are calculated from inter- *
20 * polation grids covering the regions *
21 * Q^2/GeV^2 between 0.8 and 1.E6 ( 1.E4 for F2 ) *
22 * x between 1.E-9 and 1. *
23 * Any call outside these regions stops the program execution. *
24 * *
25 * At Q^2 = MZ^2, alpha_s reads 0.114 (0.125) in NLO (LO); the *
26 * heavy quark thresholds, QH^2 = mh^2, in the beta function are *
27 * mc = 1.4 GeV, mb = 4.5 GeV, mt = 175 GeV. *
28 * Note that the NLO alpha_s running is different from GRV(94). *
29 * *
30 * Questions, comments etc to: avogt@physik.uni-wuerzburg.de *
31 * *
32 * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
33 *
34 *
35 *
36 *
37  SUBROUTINE grv98pa (ISET, X, Q2, UV, DV, US, DS, SS, GL)
38 *********************************************************************
39 * *
40 * THE PARTON ROUTINE. *
41 * __ *
42 * INPUT: ISET = 1 (LO), 2 (NLO, MS), or 3 (NLO, DIS) *
43 * X = Bjorken-x (between 1.E-9 and 1.) *
44 * Q2 = scale in GeV**2 (between 0.8 and 1.E6) *
45 * *
46 * OUTPUT: UV = u - u(bar), DV = d - d(bar), US = u(bar), *
47 * DS = d(bar), SS = s = s(bar), GL = gluon. *
48 * Always x times the distribution is returned. *
49 * *
50 * COMMON: The main program or the calling routine has to have *
51 * a common block COMMON / INTINIP / IINIP , and the *
52 * integer variable IINIP has always to be zero when *
53 * GRV98PA is called for the first time or when ISET *
54 * has been changed. *
55 * *
56 * GRIDS: 1. grv98lo.grid, 2. grv98nlm.grid, 3. grv98nld.grid, *
57 * (1+1809 lines with 6 columns, 4 significant figures) *
58 * *
59 *******************************************************i*************
60 *
61  IMPLICIT DOUBLE PRECISION (a-h, o-z)
62  parameter(npart=6, nx=68, nq=27, narg=2)
63  dimension xuvf(nx,nq), xdvf(nx,nq), xdef(nx,nq), xudf(nx,nq),
64  1 xsf(nx,nq), xgf(nx,nq), parton(npart,nq,nx-1),
65  2 qs(nq), xb(nx), xt(narg), na(narg), arrf(nx+nq)
66  CHARACTER*80 line
67  COMMON / intinip / iinip
68  SAVE xuvf, xdvf, xdef, xudf, xsf, xgf, na, arrf
69  include "pepadm.inc"
70 !bs>
71 *
72 *...BJORKEN-X AND Q**2 VALUES OF THE GRID :
73  DATA qs / 0.8e0,
74  1 1.0e0, 1.3e0, 1.8e0, 2.7e0, 4.0e0, 6.4e0,
75  2 1.0e1, 1.6e1, 2.5e1, 4.0e1, 6.4e1,
76  3 1.0e2, 1.8e2, 3.2e2, 5.7e2,
77  4 1.0e3, 1.8e3, 3.2e3, 5.7e3,
78  5 1.0e4, 2.2e4, 4.6e4,
79  6 1.0e5, 2.2e5, 4.6e5,
80  7 1.e6 /
81  DATA xb / 1.0e-9, 1.8e-9, 3.2e-9, 5.7e-9,
82  a 1.0e-8, 1.8e-8, 3.2e-8, 5.7e-8,
83  b 1.0e-7, 1.8e-7, 3.2e-7, 5.7e-7,
84  c 1.0e-6, 1.4e-6, 2.0e-6, 3.0e-6, 4.5e-6, 6.7e-6,
85  1 1.0e-5, 1.4e-5, 2.0e-5, 3.0e-5, 4.5e-5, 6.7e-5,
86  2 1.0e-4, 1.4e-4, 2.0e-4, 3.0e-4, 4.5e-4, 6.7e-4,
87  3 1.0e-3, 1.4e-3, 2.0e-3, 3.0e-3, 4.5e-3, 6.7e-3,
88  4 1.0e-2, 1.4e-2, 2.0e-2, 3.0e-2, 4.5e-2, 0.06, 0.08,
89  5 0.1, 0.125, 0.15, 0.175, 0.2, 0.225, 0.25, 0.275,
90  6 0.3, 0.325, 0.35, 0.375, 0.4, 0.45, 0.5, 0.55,
91  7 0.6, 0.65, 0.7, 0.75, 0.8, 0.85, 0.9, 0.95, 1. /
92 *
93 *...CHECK OF X AND Q2 VALUES :
94  IF ( (x.LT.0.99d-9) .OR. (x.GT.1.d0) ) THEN
95  WRITE(6,91)
96  91 FORMAT (2x,'PARTON INTERPOLATION: X OUT OF RANGE')
97  stop
98  ENDIF
99  IF ( (q2.LT.0.799) .OR. (q2.GT.1.01e6) ) THEN
100  WRITE(6,92)
101  92 FORMAT (2x,'PARTON INTERPOLATION: Q2 OUT OF RANGE')
102  stop
103  ENDIF
104  IF (iinip .NE. 0) goto 16
105 *
106 *...INITIALIZATION, IF REQUIRED :
107 *
108 * SELECTION AND READING OF THE GRID :
109 * (COMMENT: FIRST NUMBER IN THE FIRST LINE OF THE GRID)
110  write(*,*) "CUNPOL=", cunpol
111  IF (iset .EQ. 1) THEN
112  OPEN (iplst(2),file=cunpol,status='old') ! 7.332E-05
113  ELSE IF (iset .EQ. 2) THEN
114  OPEN (iplst(2),file=cunpol,status='old') ! 1.015E-04
115  ELSE IF (iset .EQ. 3) THEN
116  OPEN (iplst(2),file=cunpol,status='old') ! 1.238E-04
117  ELSE
118  WRITE(6,93)
119  93 FORMAT (2x,'NO OR INVALID PARTON SET CHOICE')
120  stop
121  END IF
122  iinip = 1
123  READ(iplst(2),89) line
124  89 FORMAT(a80)
125  DO 15 m = 1, nx-1
126  DO 15 n = 1, nq
127  READ(iplst(2),90) parton(1,n,m), parton(2,n,m), parton(3,n,m),
128  1 parton(4,n,m), parton(5,n,m), parton(6,n,m)
129  90 FORMAT (6(1pe10.3))
130  15 CONTINUE
131  CLOSE(iplst(2))
132 *
133 *....ARRAYS FOR THE INTERPOLATION SUBROUTINE :
134  DO 10 iq = 1, nq
135  DO 20 ix = 1, nx-1
136  xb0v = xb(ix)**0.5
137  xb0s = xb(ix)**(-0.2)
138  xb1 = 1.-xb(ix)
139  xuvf(ix,iq) = parton(1,iq,ix) / (xb1**3 * xb0v)
140  xdvf(ix,iq) = parton(2,iq,ix) / (xb1**4 * xb0v)
141  xdef(ix,iq) = parton(3,iq,ix) / (xb1**7 * xb0v)
142  xudf(ix,iq) = parton(4,iq,ix) / (xb1**7 * xb0s)
143  xsf(ix,iq) = parton(5,iq,ix) / (xb1**7 * xb0s)
144  xgf(ix,iq) = parton(6,iq,ix) / (xb1**5 * xb0s)
145  20 CONTINUE
146  xuvf(nx,iq) = 0.e0
147  xdvf(nx,iq) = 0.e0
148  xdef(nx,iq) = 0.e0
149  xudf(nx,iq) = 0.e0
150  xsf(nx,iq) = 0.e0
151  xgf(nx,iq) = 0.e0
152  10 CONTINUE
153  na(1) = nx
154  na(2) = nq
155  DO 30 ix = 1, nx
156  arrf(ix) = dlog(xb(ix))
157  30 CONTINUE
158  DO 40 iq = 1, nq
159  arrf(nx+iq) = dlog(qs(iq))
160  40 CONTINUE
161 *
162 *...CONTINUATION, IF INITIALIZATION WAS DONE PREVIOUSLY.
163 *
164  16 CONTINUE
165 *
166 *...INTERPOLATION :
167  xt(1) = dlog(x)
168  xt(2) = dlog(q2)
169  x1 = 1.- x
170  xv = x**0.5
171  xs = x**(-0.2)
172  uv = bsfint(narg,xt,na,arrf,xuvf) * x1**3 * xv
173  dv = bsfint(narg,xt,na,arrf,xdvf) * x1**4 * xv
174  de = bsfint(narg,xt,na,arrf,xdef) * x1**7 * xv
175  ud = bsfint(narg,xt,na,arrf,xudf) * x1**7 * xs
176  us = 0.5 * (ud - de)
177  ds = 0.5 * (ud + de)
178  ss = bsfint(narg,xt,na,arrf,xsf) * x1**7 * xs
179  gl = bsfint(narg,xt,na,arrf,xgf) * x1**5 * xs
180 *
181  60 RETURN
182  END