COVERITY reports many FORWARD_NULL, corrected; AliEMCALv2: initialization of fHits...
[u/mrichter/AliRoot.git] / LHAPDF / lhapdf5.3.1 / wrapgrv.f
CommitLineData
4e9e3152 1 subroutine GRVevolve(xin,qin,pdf)
2 implicit real*8 (a-h,o-z)
3 include 'parmsetup.inc'
4 PARAMETER(ngrid=2)
5 PARAMETER (NPART=6, NX=68, NQ=27, NARG=2)
6 character*16 name(nmxset)
7 integer nmem(nmxset),ndef(nmxset),mmem
8 common/NAME/name,nmem,ndef,mmem
9 DIMENSION XXUVF(0:ngrid,NX,NQ), XXDVF(0:ngrid,NX,NQ),
10 + XXDEF(0:ngrid,NX,NQ), XXUDF(0:ngrid,NX,NQ),
11 1 XXSF(0:ngrid,NX,NQ), XXGF(0:ngrid,NX,NQ),
12 + XUVF(NX,NQ), XDVF(NX,NQ),
13 + XDEF(NX,NQ), XUDF(NX,NQ),
14 1 XSF(NX,NQ), XGF(NX,NQ),
15 + PARTON (NPART,NQ,NX-1),
16 2 QS(NQ), XB(NX), XT(NARG), NA(NARG), ARRF(NX+NQ)
17 CHARACTER*80 LINE
18 dimension pdf(-6:6)
19 save
20 x=xin
21 q2=qin*qin
22*...CHECK OF X AND Q2 VALUES :
23 IF ( (X.LT.0.99D-9) .OR. (X.GT.1.D0) ) THEN
24 WRITE(6,91)
25 91 FORMAT (2X,'PARTON INTERPOLATION: X OUT OF RANGE')
26 STOP
27 ENDIF
28 IF ( (Q2.LT.0.799) .OR. (Q2.GT.1.01E6) ) THEN
29 WRITE(6,92)
30 92 FORMAT (2X,'PARTON INTERPOLATION: Q2 OUT OF RANGE')
31 STOP
32 ENDIF
33*...INTERPOLATION :
34 DO IQ=1,NQ
35 DO IX=1,NX
36 xuvf(ix,iq)=xxuvf(imem,ix,iq)
37 xdvf(ix,iq)=xxdvf(imem,ix,iq)
38 xdef(ix,iq)=xxdef(imem,ix,iq)
39 xudf(ix,iq)=xxudf(imem,ix,iq)
40 xsf(ix,iq)=xxsf(imem,ix,iq)
41 xgf(ix,iq)=xxgf(imem,ix,iq)
42 enddo
43 enddo
44 XT(1) = DLOG(X)
45 XT(2) = DLOG(Q2)
46 X1 = 1.- X
47 XV = X**0.5
48 XS = X**(-0.2)
49 UV = FINT(NARG,XT,NA,ARRF,XUVF) * X1**3 * XV
50 DV = FINT(NARG,XT,NA,ARRF,XDVF) * X1**4 * XV
51 DE = FINT(NARG,XT,NA,ARRF,XDEF) * X1**7 * XV
52 UD = FINT(NARG,XT,NA,ARRF,XUDF) * X1**7 * XS
53 US = 0.5 * (UD - DE)
54 DS = 0.5 * (UD + DE)
55 SS = FINT(NARG,XT,NA,ARRF,XSF) * X1**7 * XS
56 GL = FINT(NARG,XT,NA,ARRF,XGF) * X1**5 * XS
57*
58
59 pdf(-6) = 0.0d0
60 pdf(6) = 0.0d0
61 pdf(-5) = 0.0d0
62 pdf(5) = 0.0d0
63 pdf(-4) = 0.0d0
64 pdf(4) = 0.0d0
65 pdf(-3) = ss
66 pdf(3) = ss
67 pdf(-2) = us
68 pdf(2) = uv+us
69 pdf(-1) = ds
70 pdf(1) = dv+ds
71 pdf(0) = gl
72 return
73c
74ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
75 entry GRVread(nset)
76c
77c print *,'calling grvread'
78 read(1,*)nmem(nset),ndef(nset)
79
80 read(1,93)xb
81c print *,xb
82 read(1,93)qs
83c print *,qs
84 93 format(8e8.2)
85c
86 do ng=0,nmem(nset)
87c
88 READ(1,89) LINE
89 89 FORMAT(A80)
90c print *,line
91 DO 15 M = 1, NX-1
92 DO 15 N = 1, NQ
93 READ(1,90) PARTON(1,N,M), PARTON(2,N,M), PARTON(3,N,M),
94 1 PARTON(4,N,M), PARTON(5,N,M), PARTON(6,N,M)
95 90 FORMAT (6(1PE10.3))
96 15 CONTINUE
97*
98*....ARRAYS FOR THE INTERPOLATION SUBROUTINE :
99 DO 10 IQ = 1, NQ
100 DO 20 IX = 1, NX-1
101 XB0V = XB(IX)**0.5
102 XB0S = XB(IX)**(-0.2)
103 XB1 = 1.-XB(IX)
104 xXUVF(ng,IX,IQ) = PARTON(1,IQ,IX) / (XB1**3 * XB0V)
105 xXDVF(ng,IX,IQ) = PARTON(2,IQ,IX) / (XB1**4 * XB0V)
106 xXDEF(ng,IX,IQ) = PARTON(3,IQ,IX) / (XB1**7 * XB0V)
107 xXUDF(ng,IX,IQ) = PARTON(4,IQ,IX) / (XB1**7 * XB0S)
108 xXSF(ng,IX,IQ) = PARTON(5,IQ,IX) / (XB1**7 * XB0S)
109 xXGF(ng,IX,IQ) = PARTON(6,IQ,IX) / (XB1**5 * XB0S)
110 20 CONTINUE
111 xXUVF(ng,NX,IQ) = 0.E0
112 xXDVF(ng,NX,IQ) = 0.E0
113 xXDEF(ng,NX,IQ) = 0.E0
114 xXUDF(ng,NX,IQ) = 0.E0
115 xXSF(ng,NX,IQ) = 0.E0
116 xXGF(ng,NX,IQ) = 0.E0
117 10 CONTINUE
118 NA(1) = NX
119 NA(2) = NQ
120 DO 30 IX = 1, NX
121 ARRF(IX) = DLOG(XB(IX))
122 30 CONTINUE
123 DO 40 IQ = 1, NQ
124 ARRF(NX+IQ) = DLOG(QS(IQ))
125 40 CONTINUE
126
127 enddo
128
129 return
130c
131ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
132 entry GRValfa(alfas,qalfa)
133 call getnset(iset)
134 q2alfa = qalfa*qalfa
135 call GetOrderAsM(iset,iord)
136 nord=iord+1
137 alfas=grvals(q2alfa,nord)
138 alfas = 4.0d0*3.14159d0*alfas
139 return
140c
141ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
142 entry GRVinit(Eorder,Q2fit)
143 return
144c
145ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
146 entry GRVpdf(mem)
147 imem = mem
148 return
149c
150 1000 format(5e13.5)
151 end
152c
153 FUNCTION FINT(NARG,ARG,NENT,ENT,TABLE)
154*********************************************************************
155* *
156* THE INTERPOLATION ROUTINE (CERN LIBRARY ROUTINE E104) *
157* *
158*********************************************************************
159 IMPLICIT DOUBLE PRECISION (A-H, O-Z)
160 DIMENSION ARG(5),NENT(5),ENT(10),TABLE(10)
161 DIMENSION D(5),NCOMB(5),IENT(5)
162 KD=1
163 M=1
164 JA=1
165 DO 5 I=1,NARG
166 NCOMB(I)=1
167 JB=JA-1+NENT(I)
168 DO 2 J=JA,JB
169 IF (ARG(I).LE.ENT(J)) GO TO 3
170 2 CONTINUE
171 J=JB
172 3 IF (J.NE.JA) GO TO 4
173 J=J+1
174 4 JR=J-1
175 D(I)=(ENT(J)-ARG(I))/(ENT(J)-ENT(JR))
176 IENT(I)=J-JA
177 KD=KD+IENT(I)*M
178 M=M*NENT(I)
179 5 JA=JB+1
180 FINT=0.
181 10 FAC=1.
182 IADR=KD
183 IFADR=1
184 DO 15 I=1,NARG
185 IF (NCOMB(I).EQ.0) GO TO 12
186 FAC=FAC*(1.-D(I))
187 GO TO 15
188 12 FAC=FAC*D(I)
189 IADR=IADR-IFADR
190 15 IFADR=IFADR*NENT(I)
191 FINT=FINT+FAC*TABLE(IADR)
192 IL=NARG
193 40 IF (NCOMB(IL).EQ.0) GO TO 80
194 NCOMB(IL)=0
195 IF (IL.EQ.NARG) GO TO 10
196 IL=IL+1
197 DO 50 K=IL,NARG
198 50 NCOMB(K)=1
199 GO TO 10
200 80 IL=IL-1
201 IF(IL.NE.0) GO TO 40
202 RETURN
203 END
204
205c FUNCTION ALPHAS (Q2, NAORD)
206 FUNCTION grvals (Q2, NAORD)
207*********************************************************************
208* *
209* THE ALPHA_S ROUTINE. *
210* *
211* INPUT : Q2 = scale in GeV**2 (not too low, of course); *
212* NAORD = 1 (LO), 2 (NLO). *
213* *
214* OUTPUT: alphas_s/(4 pi) for use with the GRV(98) partons. *
215* *
216*******************************************************i*************
217*
218 IMPLICIT DOUBLE PRECISION (A - Z)
219 INTEGER NF, K, I, NAORD
220 DIMENSION LAMBDAL (3:6), LAMBDAN (3:6), Q2THR (3)
221*
222*...HEAVY QUARK THRESHOLDS AND LAMBDA VALUES :
223 DATA Q2THR / 1.960, 20.25, 30625. /
224 DATA LAMBDAL / 0.2041, 0.1750, 0.1320, 0.0665 /
225 DATA LAMBDAN / 0.2994, 0.2460, 0.1677, 0.0678 /
226*
227*...DETERMINATION OF THE APPROPRIATE NUMBER OF FLAVOURS :
228 NF = 3
229 DO 10 K = 1, 3
230 IF (Q2 .GT. Q2THR (K)) THEN
231 NF = NF + 1
232 ELSE
233 GO TO 20
234 END IF
235 10 CONTINUE
236*
237*...LO ALPHA_S AND BETA FUNCTION FOR NLO CALCULATION :
238 20 B0 = 11.- 2./3.* NF
239 B1 = 102.- 38./3.* NF
240 B10 = B1 / (B0*B0)
241 IF (NAORD .EQ. 1) THEN
242 LAM2 = LAMBDAL (NF) * LAMBDAL (NF)
243 ALP = 1./(B0 * DLOG (Q2/LAM2))
244 GO TO 1
245 ELSE IF (NAORD .EQ. 2) then
246 LAM2 = LAMBDAN (NF) * LAMBDAN (NF)
247 B1 = 102.- 38./3.* NF
248 B10 = B1 / (B0*B0)
249 ELSE
250 WRITE (6,91)
251 91 FORMAT ('INVALID CHOICE FOR ORDER IN ALPHA_S')
252 STOP
253 END IF
254*
255*...START VALUE FOR NLO ITERATION :
256 LQ2 = DLOG (Q2 / LAM2)
257 ALP = 1./(B0*LQ2) * (1.- B10*DLOG(LQ2)/LQ2)
258*
259*...EXACT NLO VALUE, FOUND VIA NEWTON PROCEDURE :
260 DO 2 I = 1, 6
261 XL = DLOG (1./(B0*ALP) + B10)
262 XLP = DLOG (1./(B0*ALP*1.01) + B10)
263 XLM = DLOG (1./(B0*ALP*0.99) + B10)
264 Y = LQ2 - 1./ (B0*ALP) + B10 * XL
265 Y1 = (- 1./ (B0*ALP*1.01) + B10 * XLP
266 1 + 1./ (B0*ALP*0.99) - B10 * XLP) / (0.02D0*ALP)
267 ALP = ALP - Y/Y1
268 2 CONTINUE
269*
270*...OUTPUT :
271c 1 ALPHAS = ALP
272 1 grvals = ALP
273 RETURN
274 END
275
276
277