VBF: libstructf_14.f

File libstructf_14.f, 23.0 KB (added by trac, 7 years ago)
Line 
1C------------------------------------------------------
2C------------------------------------------------------
3C------------------------------------------------------
4c------REMOVE NF FROM THE SINGLET COEFFICIENT FUNCTIONS
5C------------------------------------------------------
6C------------------------------------------------------
7C------------------------------------------------------
8
9
10c--------couplings----second argument 1->W-, 2->w+, 3->Z, 4->gamma,
11c-----                                  5->gammaZ
12      subroutine SetCouplings()
13      implicit none
14      integer nf,i,j
15      double precision couplings(-6:6,5), couplings3(-6:6,5),sthw
16      common/coupl/couplings
17      common/coupl3/couplings3
18      common/nflav/nf
19      sthw=0.2315d0
20      do i=-6,6
21       do j=1,5
22       couplings(i,j)=0d0
23       couplings3(i,j)=0d0
24       end do
25      end do
26c---------- W- ------------------     
27      do i=1,3
28      couplings(2*i,1)=2d0
29      couplings3(2*i,1)=2d0
30      couplings(-2*i+1,1)=2d0
31      couplings3(-2*i+1,1)=-2d0
32      end do
33c--------- W+ ------------------     
34      do i=-6,6
35      couplings(i,2)=couplings(-i,1)
36      couplings3(i,2)=-couplings3(-i,1)
37      end do
38c---------gamma-----------------
39      do i=1,3
40      end do
41c--------Z0--------------------
42      do i=1,3
43      couplings(2*i,3)=2d0*(1d0/4d0+(1d0/2d0-4d0/3d0*sthw)**2)
44      couplings(-2*i,3)=2d0*(1d0/4d0+(1d0/2d0-4d0/3d0*sthw)**2)
45      couplings(2*i-1,3)=2d0*(1d0/4d0+(1d0/2d0-2d0/3d0*sthw)**2)
46      couplings(-2*i+1,3)=2d0*(1d0/4d0+(1d0/2d0-2d0/3d0*sthw)**2)
47
48      couplings3(2*i,3)=(1d0-8d0/3d0*sthw)
49      couplings3(-2*i,3)=-(1d0-8d0/3d0*sthw)
50      couplings3(2*i-1,3)=(1d0-4d0/3d0*sthw)
51      couplings3(-2*i+1,3)=-(1d0-4d0/3d0*sthw)
52      end do
53      do i=1,nf
54       do j=1,5
55       couplings(0,j)=couplings(0,j)+couplings(-i,j)+couplings(i,j)
56c------the coulplings3 for the gluon is actually not needed--------
57       couplings3(0,j)=couplings3(0,j)-couplings3(-i,j)+couplings3(i,j)
58       end do
59      end do
60      return
61      end
62
63ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
64c--------------------------------------------------------------------------
65c  PDF gives the value of the parton density, i.e. q(x) instead of x q(x), with the coupling included 
66
67      subroutine PDF(x,q,f)
68      implicit none
69      integer i
70      double precision x,q,f(-6:6)
71      call Evolvepdf(x,q,f)
72      do i=-6,6
73       f(i)=f(i)/x
74       end do
75      return
76      end
77
78      double precision function pdfg(z)
79c-----the gluon pdf including the sum of the couplings   
80      implicit none
81      integer v
82      double precision z,x,q,f(-6:6),couplings(-6:6,5)
83      common/pdfpar/x,q
84      common/vect/v
85      common/coupl/couplings
86      call setcouplings()
87      call pdf(x/z,q,f)
88      pdfg=f(0)/z *couplings(0,v)
89      return
90      end
91
92      double precision function Sing(z)
93      implicit none
94      integer i,nf,v
95      common/nflav/nf
96      common/vect/v
97      double precision z,x,Q, quark, antiq, f(-6:6)
98      double precision couplings(-6:6,5)
99      common/coupl/couplings
100      common/pdfpar/x,q
101      call setcouplings()
102      call PDF(x/z,q,f)
103      quark=0d0
104      antiq=0d0
105      do i=1,nf
106      quark=quark+f(i)
107      antiq=antiq+f(-i)
108      end do
109      sing=(quark+antiq)/z *couplings(0,v)
110      return
111      end
112
113      double precision function NSing(z)
114      implicit none
115      integer i,nf,v,f3c
116      common/nflav/nf
117      common/vect/v
118      common/f3call/f3c
119      double precision z,x,Q, quark, antiq, f(-6:6)
120      double precision couplings(-6:6,5), couplings3(-6:6,5)
121      common/coupl/couplings
122      common/coupl3/couplings3
123      common/pdfpar/x,q
124      call setcouplings()
125      call PDF(x/z,q,f)
126      quark=0d0
127      antiq=0d0
128     
129      if (f3c.eq.0) then
130      do i=1,nf
131      quark=quark+f(i) *couplings(i,v)
132      antiq=antiq+f(-i) *couplings(-i,v)
133      end do
134      else if (f3c.eq.1) then
135      do i=1,nf
136      quark=quark+f(i) *couplings3(i,v)
137      antiq=antiq+f(-i) *couplings3(-i,v)
138      end do
139      endif
140
141      nsing=(quark+antiq)/z
142      return
143      end
144
145      double precision function SingReg(z)
146      implicit none
147      double precision sing,z
148      singreg=sing(z) - sing(1d0)
149      return
150      end
151
152      double precision function NSingReg(z)
153      implicit none
154      double precision nsing,z
155      nsingreg=nsing(z) - nsing(1d0)
156      return
157      end
158
159
160ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
161c--------------------------------------------------------------------------
162c------------------------   FL   ------------------------------------------
163
164c-----NLO----------------------------------------------------------
165 
166      double precision function CLNLOa(z)
167      implicit none
168      double precision z
169      clnloa=8d0/3d0 *z
170      return
171      end
172
173      double precision function cLNLOga(z)
174      implicit none
175      double precision z
176      clnloga=2d0*z*(1d0-z)
177      return
178      end
179     
180c-----NNLO---------------------------------------------------------
181c------------------------------------------------------------------
182c----------SINGLET-------------------------------------------------
183
184      DOUBLE PRECISION FUNCTION CLNNLOSA(Y)
185      IMPLICIT DOUBLE PRECISION (A-Z)
186      DL  = LOG (Y)
187      DL1 = LOG (1.D0-Y)
188      CLNNLOSA =  ( (15.94D0 - 5.212D0 * Y) * (1.D0-Y)**2 * DL1
189     1    + (0.421D0 + 1.520D0 * Y) * DL**2 + 28.09D0 * (1.D0-Y) * DL
190     2         - (2.370D0/Y - 19.27D0) * (1.D0-Y)**3 )
191      RETURN
192      END
193
194      DOUBLE PRECISION FUNCTION CLNNLOGA (Y)
195      IMPLICIT DOUBLE PRECISION (A-Z)
196      DL  = LOG (Y)
197      DL1 = LOG (1.D0-Y)
198      CLNNLOGA = ( (94.74D0 - 49.20D0 * Y) * (1.D0-Y) * DL1**2
199     1 + 864.8D0 * (1.D0-Y) * DL1 + 1161.D0* Y * DL * DL1
200     2         + 60.060 * Y * DL**2 + 39.66D0 * (1.D0-Y) * DL
201     3         - 5.333D0 * (1.D0/Y - 1.D0) )
202      RETURN
203      END
204
205c-------------------------------------------------------------------
206c---------NON SINGLET-----------------------------------------------
207C-----CHARGED CURRENT--------------     
208      DOUBLE PRECISION FUNCTION CLNNLONSA(Y)
209      IMPLICIT DOUBLE PRECISION (A-Z)
210      INTEGER NF
211      COMMON/NFLAV/NF
212      DL  = LOG (Y)
213      DL1 = LOG (1.D0-Y)
214      CLNNLONSA =
215     1   - 52.27D0 + 100.8D0 * Y
216     2   + (23.29D0 * Y - 0.043D0) * DL**2 - 22.21D0 * DL
217     3   + 13.30D0 * DL1**2 - 59.12D0 * DL1 - 141.7D0 * DL * DL1
218     4   + NF * 16.D0/27.D0 *
219     5   ( 6.D0* Y*DL1 - 12.D0* Y*DL - 25.D0* Y + 6.D0)
220      RETURN
221      END
222
223
224c--------------neutral current---------------
225      DOUBLE PRECISION  FUNCTION CLNNLONSA_nc (Y)
226      IMPLICIT DOUBLE PRECISION (A-Z)
227      INTEGER NF
228      COMMON/NFLAV/NF
229      DL  = LOG (Y)
230      DL1 = LOG (1.D0-Y)
231      CLNNLONSA_NC =
232     1          - 40.41D0 + 97.48D0
233     2          + (26.56D0 * Y - 0.031D0) * DL**2 - 14.85D0 * DL
234     3          + 13.62D0 * DL1**2 - 55.79D0 * DL1 - 150.5D0 * DL * DL1
235     4 + NF * 16.D0/27.D0 * ( 6.D0* Y*DL1 - 12.D0* Y*DL - 25.* Y + 6.D0)
236      RETURN
237      END
238
239
240      DOUBLE PRECISION FUNCTION CLNNLONSC (Y)
241      IMPLICIT DOUBLE PRECISION (A-Z)
242      integer v
243      common/vect/v
244      if (v.le.2) then
245      CLNNLONSC = -0.150D0
246      else
247      CLNNLONSC=-0.164d0
248      end if
249      RETURN
250      END
251
252      double precision function intLnloA(z)
253      implicit none
254      double precision z,cLnloa,nsing
255      intLnloa=cLnloA(z)*nsing(z)
256      return
257      end
258
259      double precision function intLnlogA(z)
260      implicit none
261      double precision z,cLnloga, pdfg
262      intLnloga=cLnloga(z)*pdfg(z)
263      return
264      end
265
266      double precision function IntLNNLOsa(z)
267      implicit none
268      double precision z, clnnlosa, sing
269      intlnnlosa=clnnlosa(z)*sing(z)
270      return
271      end
272
273      double precision function IntLNNLOga(z)
274      implicit none
275      double precision z, clnnloga, pdfg
276      intlnnloGa=clnnloGa(z)*pdfg(z)
277      return
278      end
279     
280      double precision function IntLNNLOnsa(z)
281      implicit none
282      double precision z, clnnlonsa,clnnlonsa_nc, nsing
283      integer v
284      common/vect/v
285      if(v.le.2) then
286      intlnnlonsa=clnnlonsa(z)*nsing(z)
287      else
288      intlnnlonsa=clnnlonsa_nc(z)*nsing(z)
289      end if
290      return
291      end
292
293
294      double precision function FL(x,q,ord,v,zz)
295      implicit none
296      double precision x,q,zz, nsing,alphaspdf,pi,eps
297      integer ord, nf,v,i,vv,f2,f1
298      double precision intlnloa, intlnloga
299      double precision intlnnlosa, intlnnloga, intlnnlonsa
300      double precision clnnlonsc
301      common/prec/eps
302      double precision xx,qq
303      common/pdfpar/xx,qq
304      common/nflav/nf
305      common/vect/vv
306      integer f3c
307      common/f3call/f3c
308      double precision z,z0
309       z=(1d0-eps-x)*zz +x
310       z0=zz*x
311
312      if(v.le.2) then
313      nf=4
314      else
315      nf=5
316      endif
317
318      f3c=0
319
320      xx=x
321      qq=q
322      vv=v
323      if(1d0-x.lt.eps) then
324      fl=0d0
325      else
326
327
328      pi=3.1415926535
329
330      if (ord.eq.1) then
331      FL=0d0
332
333      else if (ord.eq.2) then
334      fL=x*alphaspdf(q)/(4d0*pi)*2d0*(1d0-eps-x)*
335     1 (intlnloa(z)+intlnloga(z))
336     
337      else if (ord.eq.3) then
338      fl=x*(alphaspdf(q)/(4d0*pi))**2* (
339     1    (1d0-x-eps)*(intlnnlosa(z)+intlnnloga(z)
340     2                +intlnnlonsa(z))
341     3   +nsing(1d0)*cLnnlonsc(x)
342     4  )
343
344      endif
345      endif
346     
347      return
348      end
349     
350
351
352ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
353c--------------------------------------------------------------------------
354c------------------------   F1   ------------------------------------------
355cccc need to correct for NC
356
357
358      double precision function F1(x,q,ord,v,zz)
359      implicit none
360      double precision x,q,zz,nsing,alphaspdf,pi,eps
361      integer ord, nf,v,i,vv
362      common/prec/eps
363      double precision xx,qq
364      common/pdfpar/xx,qq
365      common/nflav/nf
366      common/vect/vv
367      double precision f2,fl
368      integer f3c
369      common/f3call/f3c
370      double precision z, z0
371      z=(1d0-eps-x)*zz +x
372      z0=zz*x
373      f3c=0
374
375      if(v.le.2) then
376      nf=4
377      else
378      nf=5
379      endif
380
381
382      xx=x
383      qq=q
384      vv=v
385      if(1d0-x.lt.eps) then
386      f1=0d0
387      else
388
389      if (ord.eq.1) then
390      F1=nsing(1d0)/2d0
391
392      else if (ord.eq.2) then
393      f1=(f2(x,q,ord,v,zz)-fl(x,q,ord,v,zz))/ (2d0*x)
394
395      else if (ord.eq.3) then
396      f1=(f2(x,q,ord,v,zz)-fl(x,q,ord,v,zz))/ (2d0*x)
397      endif
398
399      endif
400      return
401      end
402 
403
404ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
405c--------------------------------------------------------------------------
406c------------------------   F2   ------------------------------------------
407
408      double precision function C2NLOa(z)
409      implicit none
410      double precision z
411      C2NLOa= 4d0/3d0*( 3d0+2d0*z-(1d0+z)*log(1-z)
412     1   -(1d0+z**2)/(1d0-z)*log(z))
413      return
414      end
415
416      double precision function C2NLOb(z)
417      implicit none
418      double precision z
419      C2NLOb= 4d0/3d0*(2d0*log(1d0-z)/(1d0-z)-3d0/2d0 /(1d0-z) )
420      return
421      end
422
423      double precision function C2NLOc(z)
424      implicit none
425      double precision z,pi
426      pi=3.1415926535
427      C2NLOc= -4d0/3d0*(pi**2/3d0+9d0/2d0 )
428      return
429      end
430
431      double precision function C2NLOg(z)
432      implicit none
433      double precision z
434c     C2NLOg=1d0/2d0*z*( (z**2+(1-z)**2)*log((1-z)/z)-1d0+8*z*(1-z) )
435      C2NLOg=1d0/2d0*( (z**2+(1-z)**2)*log((1-z)/z)-1d0+8*z*(1-z) )
436      return
437      end
438
439
440      double precision FUNCTION C2NNLOSA(Y)
441      IMPLICIT double precision (A-Z)
442      DL  = LOG (Y)
443      DL1 = LOG (1.d0-Y)
444      C2NNLOSA =  ( 5.290d0 * (1.d0/Y-1.d0) + 4.310d0 * DL**3   
445     1   - 2.086d0 * DL**2 + 39.78d0 * DL - 0.101d0 * (1.d0-Y) * DL1**3
446     2   - (24.75d0 - 13.80d0 * Y) * DL**2 * DL1 + 30.23d0 * DL * DL1 )
447      RETURN
448      END
449
450
451      double precision FUNCTION C2NNLOGA (Y)
452      IMPLICIT double precision (A-Z)
453      DL  = LOG (Y)
454      DL1 = LOG (1.D0-Y)
455      C2NNLOGA =
456     1  ( 1.d0/Y * (11.90d0 + 1494.d0* DL1) + 5.319d0 * DL**3 
457     1       - 59.48d0 * DL**2 - 284.8d0 * DL + 392.4d0 - 1483.d0* DL1
458     2     + (6.445d0 + 209.4d0 * (1.d0-Y)) * DL1**3 - 24.00d0 * DL1**2
459     3        - 724.1d0 * DL**2 * DL1 - 871.8d0 * DL * DL1**2 )
460      RETURN
461      END
462
463      double precision FUNCTION C2NNLOGC (Y)
464      IMPLICIT double precision (A-Z)
465      C2NNLOGC=-0.28d0
466      RETURN
467      END
468c------------------NON SINGLET----------------------
469c---------------charged current-----------------------
470      double precision function C2NNLONSA(Y)
471      IMPLICIT double precision (A-Z)
472      INTEGER NF
473      COMMON/NFLAV/NF
474      DL  = LOG (Y)
475      DL1 = LOG (1.D0-Y)
476      C2NNLONSA = - 84.18d0 - 1010.d0* Y
477     2 -3.748d0 * DL**3 - 19.56d0 * DL**2 - 1.235d0 * DL
478     3 - 17.19d0 * DL1**3 + 71.08d0 * DL1**2 - 663.0d0 * DL1
479     4 - 192.4d0 * DL * DL1**2 + 80.41d0  * DL**2 * DL1
480     5 + NF * ( - 5.691d0 - 37.91d0 *Y
481     6 + 2.244d0 * DL**2 + 5.770d0 * DL
482     7 - 1.707d0* DL1**2  + 22.95d0 * DL1
483     8 + 3.036d0 * DL**2 * DL1 + 17.97d0 * DL * DL1 )     
484      RETURN
485      END
486
487      DOUBLE PRECISION FUNCTION C2NNLONSB (Y)
488      IMPLICIT DOUBLE PRECISION (A-Z)
489      INTEGER NF
490      COMMON/NFLAV/NF
491      DL1 = LOG (1.D0-Y)
492      DM  = 1./(1.D0-Y)
493      C2NNLONSB =
494     1  + 14.2222d0 * DL1**3 - 61.3333d0 * DL1**2- 31.105d0 * DL1
495     2  + 188.64d0
496     3  + NF * ( 1.77778d0 * DL1**2 - 8.5926d0 *  DL1 + 6.3489 )
497      C2NNLONSB = DM * C2NNLONSB
498      RETURN
499      END
500
501      DOUBLE PRECISION FUNCTION C2NNLONSC (Y)
502      IMPLICIT DOUBLE PRECISION (A-Z)
503      INTEGER NF
504      COMMON/NFLAV/NF
505      DL1 = LOG (1.D0-Y)
506      C2NNLONSC =
507     1 + 3.55555D0 * DL1**4 - 20.4444D0 * DL1**3 - 15.5525D0 * DL1**2
508     2 + 188.64D0 * DL1 - 338.531D0 + 0.537D0
509     3 + NF * (0.592593D0 * DL1**3 - 4.2963D0 * DL1**2
510     4 + 6.3489D0 * DL1 + 46.844D0 - 0.0035D0)
511      RETURN
512      END
513
514c------------------------ neutral current
515      DOUBLE PRECISION  FUNCTION C2NNLONSA_NC (Y)
516      IMPLICIT DOUBLE PRECISION (A-Z)
517      INTEGER NF
518      COMMON/NFLAV/NF
519      DL  = LOG (Y)
520      DL1 = LOG (1.D0-Y)
521      C2NNLONSA_NC =
522     1          - 69.59D0 - 1008.D0* Y
523     2          - 2.835D0 * DL**3 - 17.08D0 * DL**2 + 5.986D0 * DL
524     3          - 17.19D0 * DL1**3 + 71.08D0 * DL1**2 - 660.7D0 * DL1
525     4          - 174.8D0 * DL * DL1**2 + 95.09D0 * DL**2 * DL1
526     5        + NF * ( - 5.691D0 - 37.91D0 * Y
527     6          + 2.244D0 * DL**2 + 5.770D0 * DL
528     7          - 1.707D0 * DL1**2  + 22.95D0 * DL1
529     8          + 3.036D0 * DL**2 * DL1 + 17.97D0 * DL * DL1 )     
530      RETURN
531      END
532
533
534      DOUBLE PRECISION FUNCTION C2NNLONSC_NC (Y)
535      IMPLICIT REAL*8 (A-Z)
536      INTEGER NF
537      COMMON/NFLAV/NF
538      DL1 = LOG (1.D0-Y)
539      C2NNLONSC_NC =
540     1 + 3.55555D0 * DL1**4 - 20.4444D0 * DL1**3 - 15.5525D0 * DL1**2
541     2 + 188.64D0 * DL1 - 338.531D0 + 0.485D0
542     3 + NF * (0.592593D0 * DL1**3 - 4.2963D0 * DL1**2
543     4 + 6.3489D0 * DL1 + 46.844D0 - 0.0035D0)
544      RETURN
545      END
546
547
548c--------------------------------
549c--------------------------------
550c--------------------------------
551      double precision function Int2NLOa(z)
552      implicit none
553      doubLe precision z,c2nloa,nsing
554      int2nloa=c2nloa(z)*nsing(z)
555      return
556      end
557
558      double precision function Int2NLOb(z)
559      implicit none
560      doubLe precision z,c2nlob,nsingreg
561      int2nlob=c2nlob(z)*nsingreg(z)
562      return
563      end
564
565      double precision function Int2NLOg(z)
566      implicit none
567      double precision z,c2nlog,pdfg
568      int2nlog=c2nlog(z)*pdfg(z)
569      return
570      end
571
572      double precision function Int2NNLOsA(z)
573      implicit none
574      double precision z,sing, c2nnlosa
575      int2nnlosa=sing(z)*c2nnlosa(z)
576      return
577      end
578
579      double precision function Int2NNLOgA(z)
580      implicit none
581      double precision z,pdfg, c2nnloga
582      int2nnloga=pdfg(z)*c2nnloga(z)
583      return
584      end
585
586      double precision function Int2NNLOnsA(z)
587      implicit none
588      double precision z,nsing, c2nnlonsa, c2nnlonsa_nc
589      integer v
590      if (v.le.2) then
591      int2nnlonsa=nsing(z)*c2nnlonsa(z)
592      else
593      int2nnlonsa=nsing(z)*c2nnlonsa_nc(z)
594      end if
595      return
596      end
597
598      double precision function Int2NNLOnsB(z)
599      implicit none
600      double precision z,nsingreg, c2nnlonsB
601      int2nnlonsb=nsingreg(z)*c2nnlonsB(z)
602      return
603      end
604
605
606cccc need to correct for NC
607      double precision function F2(x,q,ord,v,zz)
608      implicit none
609      double precision x,q,zz,nsing,alphaspdf,pi,eps,
610     1  c2nnlogc,c2nnlonsc,c2nnlonsc_nc,c2nnlonscVAL
611      integer ord, nf,v,i,vv
612      double precision int2nloa,int2nlob, c2nlob,c2nloc, int2nlog,pdfg
613      double precision int2nnlosa,int2nnloga,
614     1 int2nnlonsa, int2nnlonsb
615      common/prec/eps
616      double precision xx,qq
617      common/pdfpar/xx,qq
618      common/nflav/nf
619      common/vect/vv
620      integer f3c
621      common/f3call/f3c
622      double precision z,z0
623      z=(1d0-eps-x)*zz +x
624      z0=zz*x
625      f3c=0
626
627      if(v.le.2) then
628      nf=4
629      else
630      nf=5
631      endif
632
633      xx=x
634      qq=q
635      vv=v
636      if(1d0-x.lt.eps) then
637      f2=0d0
638      else
639
640      pi=3.1415926535
641
642      if (ord.eq.1) then
643      F2=x*nsing(1d0)
644
645      else if (ord.eq.2) then
646      f2=alphaspdf(q)/(4d0*pi)*2d0*x*((1d0-eps-x)*
647     1 (int2nloa(z)+int2nlob(z) +int2nlog(z))
648     2 -x*nsing(1d0)*c2nlob(z0)
649     3 +nsing(1d0)*c2nloc(z)
650     1 )
651     
652      else if (ord.eq.3) then
653       if (v.le.2) then
654        c2nnlonscVAL=c2nnlonsc(x)
655       else
656        c2nnlonscVAL=c2nnlonsc_nc(x)
657       end if
658      f2=(alphaspdf(q)/(4d0*pi))**2 *x*(
659     1         (1d0-eps-x)*( int2nnlosa(z)+int2nnloga(z)
660     2                      + int2nnlonsa(z)+int2nnlonsb(z))
661     3      + pdfg(1d0)*c2nnlogc(x)
662     4      + nsing(1d0)*c2nnlonscVAL
663     5      )
664       endif
665
666      endif
667      return
668      end
669   
670 
671
672cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc--------------------------------------------------------------------------
673c------------------------   F3   ------------------------------------------
674
675      double precision function C3NLOb(z)
676      implicit none
677      double precision z,c2nlob
678      c3nlob=c2nlob(z)
679      return
680      end
681
682      double precision function C3NLOa(z)
683      implicit none
684      double precision z, c2nloa
685      c3nloa=c2nloa(z)-4d0/3d0 *(1d0+z)
686      return
687      end
688
689      DOUBLE PRECISION  FUNCTION C3NNLOSUMA (Y)
690      IMPLICIT DOUBLE PRECISION (A-Z)
691      INTEGER NF
692      COMMON/NFLAV/NF
693      DL  = LOG (Y)
694      DL1 = LOG (1.D0-Y)
695      C3NNLOSUMA =
696     1          - 206.1D0 - 576.8D0 * Y
697     2          - 3.922D0 * DL**3 - 33.31D0 * DL**2 - 67.60D0 * DL
698     3          - 15.20D0 * DL1**3 + 94.61D0 * DL1**2 - 409.6D0 * DL1
699     4          - 147.9D0 * DL * DL1**2
700     5          + NF * ( - 6.337D0 - 14.97D0 * Y
701     6          + 2.207D0 * DL**2 + 8.683D0 * DL
702     7          + 0.042D0 * DL1**3 - 0.808D0 * DL1**2 + 25.00D0 * DL1
703     8          + 9.684D0 * DL * DL1 )     
704      RETURN
705      END
706
707
708      DOUBLE PRECISION FUNCTION C3NNLODIFA (Y)
709      IMPLICIT DOUBLE PRECISION (A-Z)
710      INTEGER NF
711      COMMON/NFLAV/NF
712      DL  = LOG (Y)
713      DL1 = LOG (1.-Y)
714      C3NNLODIFA =
715     1          - 242.9D0 - 467.2D0 * Y
716     2          - 3.049D0 * DL**3 - 30.14D0 * DL**2 - 79.14D0 * DL
717     3          - 15.20D0 * DL1**3 + 94.61D0 * DL1**2 - 396.1D0 * DL1
718     4          - 92.43D0 * DL * DL1**2
719     5          + NF * ( - 6.337D0 - 14.97D0 * Y
720     6          + 2.207D0 * DL**2 + 8.683D0 * DL
721     7          + 0.042D0 * DL1**3 - 0.808D0 * DL1**2  + 25.00D0 * DL1
722     8          + 9.684D0 * DL * DL1 )     
723      RETURN
724      END
725
726      DOUBLE PRECISION  FUNCTION C3NNLOSUMC (Y)
727      IMPLICIT DOUBLE PRECISION (A-Z)
728      INTEGER NF
729      COMMON/NFLAV/NF
730      DL1 = LOG (1.D0-Y)
731      C3NNLOSUMC =
732     1    + 3.55555D0 * DL1**4 - 20.4444D0 * DL1**3 - 15.5525D0 * DL1**2
733     2        + 188.64D0 * DL1 - 338.531D0 - 0.104D0
734     3        + NF * (0.592593D0 * DL1**3 - 4.2963D0 * DL1**2
735     4        + 6.3489D0 * DL1 + 46.844D0 + 0.013D0)
736      RETURN
737      END
738
739
740      DOUBLE PRECISION FUNCTION C3NNLODIFC (Y)
741      IMPLICIT DOUBLE PRECISION (A-Z)
742      INTEGER NF
743      COMMON/NFLAV/NF
744      DL1  = LOG (1D0-Y)
745      C3NNLODIFC =
746     1  + 3.55555D0 * DL1**4 - 20.4444D0 * DL1**3 - 15.5525D0 * DL1**2
747     2        + 188.64D0 * DL1 - 338.531D0 - 0.104D0
748     3        + NF * (0.592593D0 * DL1**3 - 4.2963D0 * DL1**2
749     4        + 6.3489D0 * DL1 + 46.844D0 + 0.013D0)
750      RETURN
751      END
752
753      DOUBLE PRECISION FUNCTION C3NNLOB(Y)
754      IMPLICIT DOUBLE PRECISION (A-Z)
755      INTEGER NF
756      COMMON/NFLAV/NF
757      DL1 = LOG (1.-Y)
758      DM  = 1./(1.-Y)
759      C3NNLOB =
760     1   + 14.2222D0 * DL1**3 - 61.3333D0 * DL1**2 - 31.105D0 * DL1
761     2          + 188.64D0
762     3        + NF * ( 1.77778D0 * DL1**2 - 8.5926D0 * DL1 + 6.3489D0 )
763      C3NNLOB = DM * C3NNLOB     
764      RETURN
765      END
766
767      double precision function int3nloB(z)
768      implicit none
769      double precision z,c3nlob,nsingreg
770      int3nlob=c3nlob(z)*nsingreg(z)
771      return
772      end
773
774      double precision function int3nloA(z)
775      implicit none
776      double precision z,c3nloa,nsing
777      int3nloa=c3nloA(z)*nsing(z)
778      return
779      end
780
781      double precision function int3nnloB(z)
782      implicit none
783      double precision z,c3nnlob,nsingreg
784      int3nnlob=c3nnlob(z)*nsingreg(z)
785      return
786      end
787
788      double precision function int3nnloA(z)
789      implicit none
790      integer v
791      double precision z,c3nnloa,c3nnlosuma, c3nnlodifa,nsing
792      common/vect/v
793      if (v.eq.1) then
794      c3nnloa=(c3nnlosuma(z)+c3nnlodifa(z))/2d0
795      else if(v.eq.2) then
796      c3nnloa=(c3nnlosuma(z)-c3nnlodifa(z))/2d0
797      else if(v.eq.3) then
798      c3nnloa=c3nnlosuma(z)
799      end if
800      int3nnloa=c3nnloA*nsing(z)
801      return
802      end
803
804
805
806cccc need to correct for NC
807      double precision function F3(x,q,ord,v,zz)
808      implicit none
809      double precision x,q,zz,nsing,alphaspdf,pi,eps,
810     1 c3nnloc, c3nnlodifc, c3nnlosumc
811      integer ord, nf,v,i,vv
812      double precision int3nloa, int3nlob, c3nlob,c2nloc
813      double precision int3nnloa, int3nnlob
814      common/prec/eps
815      double precision xx,qq
816      common/pdfpar/xx,qq
817      common/nflav/nf
818      common/vect/vv
819      integer f3c
820      common/f3call/f3c
821      double precision z,z0
822      z=(1d0-eps-x)*zz +x
823      z0=zz*x
824      f3c=1
825
826      if(v.le.2) then
827      nf=4
828      else
829      nf=5
830      endif
831
832      xx=x
833      qq=q
834      vv=v
835      if(1d0-x.lt.eps) then
836      f3=0d0
837      else
838
839      if (v.eq.4) then
840      f3=0d0
841      else
842      pi=3.1415926535
843
844      if (ord.eq.1) then
845      F3=nsing(1d0)
846
847      else if (ord.eq.2) then
848       f3=alphaspdf(q)/(4d0*pi)*2d0*(
849     1  (1d0-eps-x)*(int3nloa(z)+int3nlob(z))
850     2  -x*nsing(1d0)*c3nlob(z0)
851     3  +nsing(1d0)*c2nloc(z)
852     4    )
853
854      else if (ord.eq.3) then
855       if (v.eq.1) then
856        c3nnloc=(c3nnlosumc(x)+c3nnlodifc(x))/2d0
857       else if(v.eq.2) then
858        c3nnloc=(c3nnlosumc(x)-c3nnlodifc(x))/2d0
859       else if(v.eq.3) then
860        c3nnloc=c3nnlosumc(x)
861       end if
862      f3=(alphaspdf(q)/(4d0*pi))**2 *(
863     1  (1d0-eps-x)*(int3nnloa(z)  +int3nnlob(z))
864     2 + nsing(1d0)*c3nnloc
865     3  )
866
867      endif
868     
869      endif
870      endif
871      return
872      end