VBF: libstructf_13.f

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