]>
Commit | Line | Data |
---|---|---|
1 | subroutine MRST2006evolve(x,Q,pdf) | |
2 | implicit none | |
3 | c implicit real*8 (a-h,o-z) | |
4 | include 'parmsetup.inc' | |
5 | character*16 name(nmxset) | |
6 | integer i,f,nhess,nx,nq,np,nqc0,nqb0,nqc,nqb,n,m,io | |
7 | double precision x,q,xmin,xmax,qsqmin,qsqmax,emc2,emb2,eps, | |
8 | & dummy,qsq,xlog,qsqlog,res,ExtrapLHAPDF,InterpLHAPDF | |
9 | double precision emc2minuseps,emb2minuseps | |
10 | double precision xsave,q2save | |
11 | double precision pdf(-6:6),alfas,Qalfa,Q2fit,eorder | |
12 | double precision upv,dnv,glu,usea,str,dsea,sbar,chm,bot | |
13 | integer iset,nflav,mem,nset,nmem(nmxset),ndef(nmxset),mmem,imem | |
14 | common/NAME/name,nmem,ndef,mmem | |
15 | parameter(nx=64,nq=48,np=9,nqc0=4,nqb0=14,nqc=nq-nqc0,nqb=nq-nqb0) | |
16 | parameter(xmin=1d-6,xmax=1d0,qsqmin=1d0,qsqmax=1d9,eps=1d-6) | |
17 | parameter(nhess=30,emc2=1.43d0**2,emb2=4.3d0**2) ! MRST 2006 NNLO | |
18 | parameter(emc2minuseps=emc2-eps,emb2minuseps=emb2-eps) | |
19 | double precision f1(nx,nq),f2(nx,nq),f3(nx,nq),f4(nx,nq), | |
20 | & f5(nx,nq),f6(nx,nq),f7(nx,nq),f8(nx,nq),fc(nx,nqc), | |
21 | & fb(nx,nqb),f9(nx,nq) | |
22 | double precision qq(nq),xx(nx),cc1(0:nhess,nx,nq,4,4), | |
23 | & cc2(0:nhess,nx,nq,4,4),cc3(0:nhess,nx,nq,4,4), | |
24 | & cc4(0:nhess,nx,nq,4,4),cc6(0:nhess,nx,nq,4,4), | |
25 | & cc8(0:nhess,nx,nq,4,4),cc9(0:nhess,nx,nq,4,4), | |
26 | & ccc(0:nhess,nx,nqc,4,4),ccb(0:nhess,nx,nqb,4,4) | |
27 | double precision xxl(nx),qql(nq),qqlc(nqc),qqlb(nqb) | |
28 | data xx/1d-6,2d-6,4d-6,6d-6,8d-6, | |
29 | & 1d-5,2d-5,4d-5,6d-5,8d-5, | |
30 | & 1d-4,2d-4,4d-4,6d-4,8d-4, | |
31 | & 1d-3,2d-3,4d-3,6d-3,8d-3, | |
32 | & 1d-2,1.4d-2,2d-2,3d-2,4d-2,6d-2,8d-2, | |
33 | & .1d0,.125d0,.15d0,.175d0,.2d0,.225d0,.25d0,.275d0, | |
34 | & .3d0,.325d0,.35d0,.375d0,.4d0,.425d0,.45d0,.475d0, | |
35 | & .5d0,.525d0,.55d0,.575d0,.6d0,.625d0,.65d0,.675d0, | |
36 | & .7d0,.725d0,.75d0,.775d0,.8d0,.825d0,.85d0,.875d0, | |
37 | & .9d0,.925d0,.95d0,.975d0,1d0/ | |
38 | data qq/1.d0, | |
39 | & 1.25d0,1.5d0,emc2minuseps,emc2,2.5d0,3.2d0,4d0,5d0,6.4d0,8d0, | |
40 | & 1d1,1.2d1,emb2minuseps,emb2,2.6d1,4d1,6.4d1,1d2, | |
41 | & 1.6d2,2.4d2,4d2,6.4d2,1d3,1.8d3,3.2d3,5.6d3,1d4, | |
42 | & 1.8d4,3.2d4,5.6d4,1d5,1.8d5,3.2d5,5.6d5,1d6, | |
43 | & 1.8d6,3.2d6,5.6d6,1d7,1.8d7,3.2d7,5.6d7,1d8, | |
44 | & 1.8d8,3.2d8,5.6d8,1d9/ | |
45 | save | |
46 | ||
47 | xsave=x | |
48 | qsq=q*q | |
49 | q2save=qsq | |
50 | ||
51 | ||
52 | if (qsq.gt.qq(nqc0).and.qsq.lt.qq(nqc0+1)) qsq = qq(nqc0) | |
53 | if (qsq.gt.qq(nqb0).and.qsq.lt.qq(nqb0+1)) qsq = qq(nqb0) | |
54 | ||
55 | xlog=log10(x) | |
56 | qsqlog=log10(qsq) | |
57 | ||
58 | call getnset(iset) | |
59 | call getnmem(iset,imem) | |
60 | ||
61 | upv = 0.d0 | |
62 | dnv = 0.d0 | |
63 | glu = 0.d0 | |
64 | usea = 0.d0 | |
65 | str = 0.d0 | |
66 | dsea = 0.d0 | |
67 | sbar = 0.d0 | |
68 | chm = 0.d0 | |
69 | bot = 0.d0 | |
70 | ||
71 | if (x.le.0.d0.or.x.gt.xmax.or.qsq.lt.qsqmin) then | |
72 | print *,"Error in GetOnePDF: x,qsq = ",x,qsq | |
73 | else if (x.lt.xmin.or.qsq.gt.qsqmax) then ! extrapolate | |
74 | print *, "Warning in GetOnePDF, extrapolating: f = ",f, | |
75 | & ", x = ",x,", q = ",q | |
76 | upv = ExtrapLHAPDF(imem,nhess,xlog,qsqlog,nx,nq,xxl,qql,cc1) | |
77 | dnv = ExtrapLHAPDF(imem,nhess,xlog,qsqlog,nx,nq,xxl,qql,cc2) | |
78 | glu = ExtrapLHAPDF(imem,nhess,xlog,qsqlog,nx,nq,xxl,qql,cc3) | |
79 | usea = ExtrapLHAPDF(imem,nhess,xlog,qsqlog,nx,nq,xxl,qql,cc4) | |
80 | str = ExtrapLHAPDF(imem,nhess,xlog,qsqlog,nx,nq,xxl,qql,cc6) | |
81 | dsea = ExtrapLHAPDF(imem,nhess,xlog,qsqlog,nx,nq,xxl,qql,cc8) | |
82 | sbar = ExtrapLHAPDF(imem,nhess,xlog,qsqlog,nx,nq,xxl,qql,cc9) | |
83 | if (qsq.ge.emc2) then ! chm | |
84 | chm = ExtrapLHAPDF(imem,nhess,xlog,qsqlog,nx,nqc,xxl,qqlc,ccc) | |
85 | endif | |
86 | if (qsq.ge.emb2) then ! bot | |
87 | bot = ExtrapLHAPDF(imem,nhess,xlog,qsqlog,nx,nqb,xxl,qqlb,ccb) | |
88 | endif | |
89 | else | |
90 | upv = InterpLHAPDF(imem,nhess,xlog,qsqlog,nx,nq,xxl,qql,cc1) | |
91 | dnv = InterpLHAPDF(imem,nhess,xlog,qsqlog,nx,nq,xxl,qql,cc2) | |
92 | glu = InterpLHAPDF(imem,nhess,xlog,qsqlog,nx,nq,xxl,qql,cc3) | |
93 | usea = InterpLHAPDF(imem,nhess,xlog,qsqlog,nx,nq,xxl,qql,cc4) | |
94 | str = InterpLHAPDF(imem,nhess,xlog,qsqlog,nx,nq,xxl,qql,cc6) | |
95 | dsea = InterpLHAPDF(imem,nhess,xlog,qsqlog,nx,nq,xxl,qql,cc8) | |
96 | sbar = InterpLHAPDF(imem,nhess,xlog,qsqlog,nx,nq,xxl,qql,cc9) | |
97 | if (qsq.ge.emc2) then ! chm | |
98 | chm = InterpLHAPDF(imem,nhess,xlog,qsqlog,nx,nqc,xxl,qqlc,ccc) | |
99 | endif | |
100 | if (qsq.ge.emb2) then ! bot | |
101 | bot = InterpLHAPDF(imem,nhess,xlog,qsqlog,nx,nqb,xxl,qqlb,ccb) | |
102 | endif | |
103 | ||
104 | end if | |
105 | ||
106 | pdf(0) = glu | |
107 | pdf(1) = dnv+dsea | |
108 | pdf(-1) = dsea | |
109 | pdf(2) = upv+usea | |
110 | pdf(-2) = usea | |
111 | pdf(3) = str | |
112 | pdf(-3) = sbar | |
113 | pdf(4) = chm | |
114 | pdf(-4) = chm | |
115 | pdf(5) = bot | |
116 | pdf(-5) = bot | |
117 | pdf(6) = 0.0d0 | |
118 | pdf(-6) = 0.0d0 | |
119 | ||
120 | x=xsave | |
121 | qsq=q2save | |
122 | return | |
123 | * | |
124 | ||
125 | C---------------------------------------------------------------------- | |
126 | ||
127 | entry MRST2006read(nset) | |
128 | read(1,*)nmem(nset),ndef(nset) | |
129 | ||
130 | ||
131 | do i = 0,nhess | |
132 | ||
133 | ||
134 | do n=1,nx-1 | |
135 | do m=1,nq | |
136 | c print *,i,n,m | |
137 | read(1,50) f1(n,m),f2(n,m),f3(n,m),f4(n,m), | |
138 | & f5(n,m),f7(n,m),f6(n,m),f8(n,m),f9(n,m) | |
139 | c write(6,50) f1(n,m),f2(n,m),f3(n,m),f4(n,m), | |
140 | c & f5(n,m),f7(n,m),f6(n,m),f8(n,m),f9(n,m) | |
141 | C-- Notation:1=upv 2=dnv 3=glu 4=usea 5=chm 6=str 7=bot 8=dsea 9=sbar | |
142 | enddo | |
143 | enddo | |
144 | c write(*,*)'PDF set ',i,' first element ',f1(1,1) | |
145 | do m=1,nq | |
146 | f1(nx,m)=0d0 | |
147 | f2(nx,m)=0d0 | |
148 | f3(nx,m)=0d0 | |
149 | f4(nx,m)=0d0 | |
150 | f5(nx,m)=0d0 | |
151 | f6(nx,m)=0d0 | |
152 | f7(nx,m)=0d0 | |
153 | f8(nx,m)=0d0 | |
154 | f9(nx,m)=0d0 | |
155 | enddo | |
156 | ||
157 | do n=1,nx | |
158 | xxl(n)=log10(xx(n)) | |
159 | enddo | |
160 | do m=1,nq | |
161 | qql(m)=log10(qq(m)) | |
162 | enddo | |
163 | ||
164 | call InitialLHAPDF(i,nhess,nx,nq,nqc0,nqb0,xxl,qql,f1,cc1) | |
165 | call InitialLHAPDF(i,nhess,nx,nq,nqc0,nqb0,xxl,qql,f2,cc2) | |
166 | call InitialLHAPDF(i,nhess,nx,nq,nqc0,nqb0,xxl,qql,f3,cc3) | |
167 | call InitialLHAPDF(i,nhess,nx,nq,nqc0,nqb0,xxl,qql,f4,cc4) | |
168 | call InitialLHAPDF(i,nhess,nx,nq,nqc0,nqb0,xxl,qql,f6,cc6) | |
169 | call InitialLHAPDF(i,nhess,nx,nq,nqc0,nqb0,xxl,qql,f8,cc8) | |
170 | call InitialLHAPDF(i,nhess,nx,nq,nqc0,nqb0,xxl,qql,f9,cc9) | |
171 | ||
172 | do m=1,nqc | |
173 | qqlc(m)=qql(m+nqc0) | |
174 | do n=1,nx | |
175 | fc(n,m)=f5(n,m+nqc0) | |
176 | enddo | |
177 | enddo | |
178 | call InitialLHAPDF(i,nhess,nx,nqc,nqc0-nqc0,nqb0-nqc0, | |
179 | & xxl,qqlc,fc,ccc) | |
180 | ||
181 | do m=1,nqb | |
182 | qqlb(m)=qql(m+nqb0) | |
183 | do n=1,nx | |
184 | fb(n,m)=f7(n,m+nqb0) | |
185 | enddo | |
186 | enddo | |
187 | call InitialLHAPDF(i,nhess,nx,nqb,nqc0-nqb0,nqb0-nqb0, | |
188 | & xxl,qqlb,fb,ccb) | |
189 | ||
190 | enddo | |
191 | ||
192 | 50 format(9e13.5) | |
193 | return | |
194 | C---------------------------------------------------------------------- | |
195 | ||
196 | entry MRST2006alfa(nflav,alfas,Qalfa) | |
197 | c print *,'MRST2006alfa',nflav,alfas,Qalfa | |
198 | call getnset(iset) | |
199 | call alphamrs(nflav,alfas,Qalfa) | |
200 | c print *,'MRST2006alfa',nflav,alfas,Qalfa | |
201 | return | |
202 | * | |
203 | C---------------------------------------------------------------------- | |
204 | entry MRST2006init(Eorder,Q2fit) | |
205 | return | |
206 | C---------------------------------------------------------------------- | |
207 | * | |
208 | entry MRST2006pdf(mem) | |
209 | call getnset(iset) | |
210 | call setnmem(iset,mem) | |
211 | return | |
212 | * | |
213 | end | |
214 | ||
215 | C---------------------------------------------------------------------- | |
216 | ||
217 | subroutine InitialLHAPDF(i,nhess,nx,my,myc0,myb0,xx,yy,ff,cc) | |
218 | implicit none | |
219 | integer nhess,i,nx,my,myc0,myb0,j,k,l,m,n | |
220 | double precision xx(nx),yy(my),ff(nx,my), | |
221 | & ff1(nx,my),ff2(nx,my), | |
222 | & ff12(nx,my),yy0(4),yy1(4),yy2(4),yy12(4),z(16), | |
223 | & cl(16),cc(0:nhess,nx,my,4,4),iwt(16,16), | |
224 | & dx,dy,polderiv,d1,d2,d1d2,xxd | |
225 | ||
226 | data iwt/1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0, | |
227 | & 0,0,0,0,0,0,0,0,1,0,0,0,0,0,0,0, | |
228 | & -3,0,0,3,0,0,0,0,-2,0,0,-1,0,0,0,0, | |
229 | & 2,0,0,-2,0,0,0,0,1,0,0,1,0,0,0,0, | |
230 | & 0,0,0,0,1,0,0,0,0,0,0,0,0,0,0,0, | |
231 | & 0,0,0,0,0,0,0,0,0,0,0,0,1,0,0,0, | |
232 | & 0,0,0,0,-3,0,0,3,0,0,0,0,-2,0,0,-1, | |
233 | & 0,0,0,0,2,0,0,-2,0,0,0,0,1,0,0,1, | |
234 | & -3,3,0,0,-2,-1,0,0,0,0,0,0,0,0,0,0, | |
235 | & 0,0,0,0,0,0,0,0,-3,3,0,0,-2,-1,0,0, | |
236 | & 9,-9,9,-9,6,3,-3,-6,6,-6,-3,3,4,2,1,2, | |
237 | & -6,6,-6,6,-4,-2,2,4,-3,3,3,-3,-2,-1,-1,-2, | |
238 | & 2,-2,0,0,1,1,0,0,0,0,0,0,0,0,0,0, | |
239 | & 0,0,0,0,0,0,0,0,2,-2,0,0,1,1,0,0, | |
240 | & -6,6,-6,6,-3,-3,3,3,-4,4,2,-2,-2,-2,-1,-1, | |
241 | & 4,-4,4,-4,2,2,-2,-2,2,-2,-2,2,1,1,1,1/ | |
242 | ||
243 | do m=1,my | |
244 | dx=xx(2)-xx(1) | |
245 | ff1(1,m)=(ff(2,m)-ff(1,m))/dx | |
246 | dx=xx(nx)-xx(nx-1) | |
247 | ff1(nx,m)=(ff(nx,m)-ff(nx-1,m))/dx | |
248 | do n=2,nx-1 | |
249 | ff1(n,m)=polderiv(xx(n-1),xx(n),xx(n+1),ff(n-1,m),ff(n,m), | |
250 | x ff(n+1,m)) | |
251 | enddo | |
252 | enddo | |
253 | ||
254 | C-- Calculate the derivatives at qsq=emc2-eps,emc2,emb2-eps,emb2 | |
255 | C-- in a similar way as at the endpoints qsqmin and qsqmax. | |
256 | do n=1,nx | |
257 | do m = 1, my | |
258 | if (m.eq.1.or.m.eq.myc0+1.or.m.eq.myb0+1) then | |
259 | dy=yy(m+1)-yy(m) | |
260 | ff2(n,m)=(ff(n,m+1)-ff(n,m))/dy | |
261 | else if (m.eq.my.or.m.eq.myc0.or.m.eq.myb0) then | |
262 | dy=yy(m)-yy(m-1) | |
263 | ff2(n,m)=(ff(n,m)-ff(n,m-1))/dy | |
264 | else | |
265 | ff2(n,m)=polderiv(yy(m-1),yy(m),yy(m+1),ff(n,m-1), | |
266 | & ff(n,m),ff(n,m+1)) | |
267 | end if | |
268 | end do | |
269 | end do | |
270 | ||
271 | do m=1,my | |
272 | dx=xx(2)-xx(1) | |
273 | ff12(1,m)=(ff2(2,m)-ff2(1,m))/dx | |
274 | dx=xx(nx)-xx(nx-1) | |
275 | ff12(nx,m)=(ff2(nx,m)-ff2(nx-1,m))/dx | |
276 | do n=2,nx-1 | |
277 | ff12(n,m)=polderiv(xx(n-1),xx(n),xx(n+1),ff2(n-1,m), | |
278 | x ff2(n,m),ff2(n+1,m)) | |
279 | enddo | |
280 | enddo | |
281 | ||
282 | do n=1,nx-1 | |
283 | do m=1,my-1 | |
284 | d1=xx(n+1)-xx(n) | |
285 | d2=yy(m+1)-yy(m) | |
286 | d1d2=d1*d2 | |
287 | ||
288 | yy0(1)=ff(n,m) | |
289 | yy0(2)=ff(n+1,m) | |
290 | yy0(3)=ff(n+1,m+1) | |
291 | yy0(4)=ff(n,m+1) | |
292 | ||
293 | yy1(1)=ff1(n,m) | |
294 | yy1(2)=ff1(n+1,m) | |
295 | yy1(3)=ff1(n+1,m+1) | |
296 | yy1(4)=ff1(n,m+1) | |
297 | ||
298 | yy2(1)=ff2(n,m) | |
299 | yy2(2)=ff2(n+1,m) | |
300 | yy2(3)=ff2(n+1,m+1) | |
301 | yy2(4)=ff2(n,m+1) | |
302 | ||
303 | yy12(1)=ff12(n,m) | |
304 | yy12(2)=ff12(n+1,m) | |
305 | yy12(3)=ff12(n+1,m+1) | |
306 | yy12(4)=ff12(n,m+1) | |
307 | ||
308 | do k=1,4 | |
309 | z(k)=yy0(k) | |
310 | z(k+4)=yy1(k)*d1 | |
311 | z(k+8)=yy2(k)*d2 | |
312 | z(k+12)=yy12(k)*d1d2 | |
313 | enddo | |
314 | ||
315 | do l=1,16 | |
316 | xxd=0.d0 | |
317 | do k=1,16 | |
318 | xxd=xxd+iwt(k,l)*z(k) | |
319 | enddo | |
320 | cl(l)=xxd | |
321 | enddo | |
322 | l=0 | |
323 | do k=1,4 | |
324 | do j=1,4 | |
325 | l=l+1 | |
326 | cc(i,n,m,k,j)=cl(l) | |
327 | enddo | |
328 | enddo | |
329 | enddo | |
330 | enddo | |
331 | return | |
332 | end | |
333 | ||
334 | C---------------------------------------------------------------------- | |
335 | ||
336 | double precision function InterpLHAPDF(i,nhess,x,y,nx,my,xx,yy, | |
337 | & cc) | |
338 | implicit none | |
339 | integer i,nx,my,nhess,locx,l,m,n | |
340 | double precision xx(nx),yy(my),cc(0:nhess,nx,my,4,4), | |
341 | & x,y,z,t,u | |
342 | ||
343 | n=locx(xx,nx,x) | |
344 | m=locx(yy,my,y) | |
345 | ||
346 | t=(x-xx(n))/(xx(n+1)-xx(n)) | |
347 | u=(y-yy(m))/(yy(m+1)-yy(m)) | |
348 | ||
349 | z=0.d0 | |
350 | do l=4,1,-1 | |
351 | z=t*z+((cc(i,n,m,l,4)*u+cc(i,n,m,l,3))*u | |
352 | . +cc(i,n,m,l,2))*u+cc(i,n,m,l,1) | |
353 | enddo | |
354 | ||
355 | InterpLHAPDF = z | |
356 | ||
357 | return | |
358 | end | |
359 | ||
360 | C---------------------------------------------------------------------- | |
361 | ||
362 | double precision function ExtrapLHAPDF(i,nhess,x,y,nx,my,xx,yy, | |
363 | & cc) | |
364 | implicit none | |
365 | integer i,nx,my,nhess,locx,n,m | |
366 | double precision xx(nx),yy(my),cc(0:nhess,nx,my,4,4), | |
367 | & x,y,z,f0,f1,z0,z1,InterpLHAPDF | |
368 | ||
369 | n=locx(xx,nx,x) ! 0: below xmin, nx: above xmax | |
370 | m=locx(yy,my,y) ! 0: below qsqmin, my: above qsqmax | |
371 | ||
372 | C-- If extrapolation in small x only: | |
373 | if (n.eq.0.and.m.gt.0.and.m.lt.my) then | |
374 | f0 = InterpLHAPDF(i,nhess,xx(1),y,nx,my,xx,yy,cc) | |
375 | f1 = InterpLHAPDF(i,nhess,xx(2),y,nx,my,xx,yy,cc) | |
376 | if (f0.gt.0.d0.and.f1.gt.0.d0) then | |
377 | z = exp(log(f0)+(log(f1)-log(f0))/(xx(2)-xx(1))*(x-xx(1))) | |
378 | else | |
379 | z = f0+(f1-f0)/(xx(2)-xx(1))*(x-xx(1)) | |
380 | end if | |
381 | C-- If extrapolation into large q only: | |
382 | else if (n.gt.0.and.m.eq.my) then | |
383 | f0 = InterpLHAPDF(i,nhess,x,yy(my),nx,my,xx,yy,cc) | |
384 | f1 = InterpLHAPDF(i,nhess,x,yy(my-1),nx,my,xx,yy,cc) | |
385 | if (f0.gt.0.d0.and.f1.gt.0.d0) then | |
386 | z = exp(log(f0)+(log(f0)-log(f1))/(yy(my)-yy(my-1))* | |
387 | & (y-yy(my))) | |
388 | else | |
389 | z = f0+(f0-f1)/(yy(my)-yy(my-1))*(y-yy(my)) | |
390 | end if | |
391 | C-- If extrapolation into large q AND small x: | |
392 | else if (n.eq.0.and.m.eq.my) then | |
393 | f0 = InterpLHAPDF(i,nhess,xx(1),yy(my),nx,my,xx,yy,cc) | |
394 | f1 = InterpLHAPDF(i,nhess,xx(1),yy(my-1),nx,my,xx,yy,cc) | |
395 | if (f0.gt.0.d0.and.f1.gt.0.d0) then | |
396 | z0 = exp(log(f0)+(log(f0)-log(f1))/(yy(my)-yy(my-1))* | |
397 | & (y-yy(my))) | |
398 | else | |
399 | z0 = f0+(f0-f1)/(yy(my)-yy(my-1))*(y-yy(my)) | |
400 | end if | |
401 | f0 = InterpLHAPDF(i,nhess,xx(2),yy(my),nx,my,xx,yy,cc) | |
402 | f1 = InterpLHAPDF(i,nhess,xx(2),yy(my-1),nx,my,xx,yy,cc) | |
403 | if (f0.gt.0.d0.and.f1.gt.0.d0) then | |
404 | z1 = exp(log(f0)+(log(f0)-log(f1))/(yy(my)-yy(my-1))* | |
405 | & (y-yy(my))) | |
406 | else | |
407 | z1 = f0+(f0-f1)/(yy(my)-yy(my-1))*(y-yy(my)) | |
408 | end if | |
409 | if (z0.gt.0.d0.and.z1.gt.0.d0) then | |
410 | z = exp(log(z0)+(log(z1)-log(z0))/(xx(2)-xx(1))*(x-xx(1))) | |
411 | else | |
412 | z = z0+(z1-z0)/(xx(2)-xx(1))*(x-xx(1)) | |
413 | end if | |
414 | else | |
415 | print *,"Error in ExtrapLHAPDF" | |
416 | stop | |
417 | end if | |
418 | ||
419 | ExtrapLHAPDF = z | |
420 | ||
421 | return | |
422 | end | |
423 | ||
424 | C---------------------------------------------------------------------- |