Warning, file /LQGENEP/lqgenep.f was not indexed
or was modified since last indexation (in which case cross-reference links may be missing, inaccurate or erroneous).
0001
0002
0003
0004
0005 subroutine LQGENEP(Nevt,flag)
0006
0007
0008
0009
0010
0011 IMPLICIT DOUBLE PRECISION(A-H, O-Z)
0012 Integer flag, itau,id,ii,ij
0013
0014
0015 double precision BEAMPAR,LQGPAR3,
0016 > ptnt,phnt,ptt,pht,pth,phh,
0017 > ptx,pty,ptz,phx,phy,phz,
0018 > ptid, phid,ppid,pxid,pyid,pzid
0019
0020 CHARACTER*256 out_file
0021 integer LQGPAR1,LQGPAR2
0022 COMMON/LQGDATA/BEAMPAR(3),LQGPAR1(10),LQGPAR2(10),LQGPAR3(20),
0023 > out_file
0024
0025
0026 double precision LQGKPAR,LQGDDX
0027 integer LQGPID
0028 COMMON/LQGEVT/LQGKPAR(3),LQGDDX(3),LQGPID(3)
0029
0030
0031
0032 INTEGER PYK,PYCHGE,PYCOMP
0033
0034
0035 PARAMETER (KSUSY1=1000000,KSUSY2=2000000,KEXCIT=4000000)
0036
0037
0038 EXTERNAL PYDATA
0039
0040
0041
0042
0043 COMMON/PYJETS/N,NPAD,K(4000,5),P(4000,5),V(4000,5)
0044
0045 COMMON/PYDAT1/MSTU(200),PARU(200),MSTJ(200),PARJ(200)
0046
0047 COMMON/PYDAT2/KCHG(500,4),PMAS(500,4),PARF(2000),VCKM(4,4)
0048
0049 COMMON/PYDAT3/MDCY(500,3),MDME(4000,2),BRAT(4000),KFDP(4000,5)
0050
0051 COMMON/PYSUBS/MSEL,MSELPD,MSUB(500),KFIN(2,-40:40),CKIN(200)
0052
0053 COMMON/PYPARS/MSTP(200),PARP(200),MSTI(200),PARI(200)
0054
0055 COMMON/PYINT1/MINT(400),VINT(400)
0056 COMMON/PYINT2/ISET(500),KFPR(500,2),COEF(500,20),ICOL(40,4,2)
0057
0058 COMMON/PYMSSM/IMSS(0:99),RMSS(0:99)
0059 SAVE /PYJETS/,/PYDAT1/,/PYDAT2/,/PYDAT3/,/PYSUBS/,/PYPARS/,
0060 >/PYINT2/,/PYMSSM/
0061
0062 DIMENSION NCHN(12),QVEC(4)
0063 DATA NCHN/12*0/
0064
0065
0066
0067 integer echar
0068 double precision ebeam,pbeam
0069 common /LQGbeam/ebeam,pbeam,echar
0070
0071
0072 character*20 parm(20)
0073 double precision pdfsf(20)
0074 common /LQGpdfC/ pdfsf
0075
0076
0077 double precision xmax,xmin,ymax,ymin,zmax,zmin,Q2min
0078 common /LQGKinC/ xmax,xmin,ymax,ymin,zmax,zmin,Q2min
0079
0080
0081 double precision Mlq,G1,G2
0082 Integer LQtype,l_o,q_i,q_j
0083 common /LQGproc/ Mlq,G1,G2,LQtype,l_o,q_i,q_j
0084
0085
0086 double precision S,Srad,x,y,z,Q2
0087 common /LQGKinV/ S,Srad,x,y,z,Q2
0088
0089
0090 double precision DXSEC(3),pvalence
0091 integer q_o,q_s,genproc,genprtyp,sch,uch,gproc(8)
0092 common /LQGout/ DXSEC,pvalence,q_o,q_s,genproc,genprtyp,
0093 >sch,uch,gproc
0094
0095
0096
0097 integer ibea
0098 Character*9 chbea(2)
0099 Character*12 chprod(2,3)
0100 data chbea /'e+ qi -> ','e- qi -> '/
0101 data chprod /' -> e+ qj',' -> e- qj',
0102 > ' -> mu+ qj',' -> mu- qj',
0103 > ' -> tau+ qj',' -> tau- qj'/
0104
0105 Character*7 LQCHA(14)
0106 DATA LQCHA /'S_0L','S_0R','~S_0R','S_1L',
0107 > 'V_1/2L','V_1/2R','~V_1/2L',
0108 > 'V_0L','V_0R','~V_0R','V_1L',
0109 > 'S_1/2L','S_1/2R','~S_1/2L'/
0110
0111
0112 Real pxtot,pytot,pztot,etot,ecmtot,
0113 > pxsum,pysum,pzsum,esum,ecmsum,
0114 > pxlow,pxhi,pylow,pyhi,pzlow,pzhi,elow,ehi,ecmlow,ecmhi
0115
0116
0117 INTEGER asciiN
0118 PARAMETER (asciiN=29)
0119 ! CHARACTER*256 out_file
0120
0121
0122
0123 Integer Nwds_HBOOK
0124 Parameter (Nwds_HBOOK=100000)
0125 Real HMEM
0126 Common /PAWC/ HMEM(Nwds_HBOOK)
0127
0128
0129
0130 ! out_file='./TestOut.txt'
0131
0132
0133
0134 If(flag.eq.0)then
0135
0136
0137 call LQGBAN
0138
0139
0140 if(LQGPAR1(3).gt.0)Call HLIMIT(Nwds_HBOOK)
0141 if(LQGPAR1(3).lt.0)Call HLIMIT(-Nwds_HBOOK)
0142
0143
0144 echar=beampar(1)
0145 Ebeam=beampar(2)
0146 Pbeam=beampar(3)
0147 S=4.d0*Ebeam*Pbeam
0148
0149
0150 MLQ=LQGPAR3(1)
0151 G1=LQGPAR3(2)
0152 G2=LQGPAR3(3)
0153 LQTYPE=LQGPAR2(1)
0154
0155
0156 l_o=LQGPAR2(4)
0157
0158
0159 q_i=LQGPAR2(2)
0160 q_j=LQGPAR2(3)
0161
0162
0163 xmin=LQGPAR3(4)
0164 xmax=LQGPAR3(5)
0165 ymin=LQGPAR3(6)
0166 ymax=LQGPAR3(7)
0167 Zmin=0.d0
0168 Zmax=1.d0
0169 Q2min=LQGPAR3(8)
0170
0171
0172 call LQGPRI1
0173
0174
0175 parm(1)='NPTYPE'
0176 parm(2)='NGROUP'
0177 parm(3)='NSET'
0178 pdfsf(1)= LQGPAR3(9)
0179 pdfsf(2)= LQGPAR3(10)
0180 pdfsf(3)= LQGPAR3(11)
0181 call PDFSET(PARM,pdfsf)
0182
0183
0184 P(1,1)=0D0
0185 P(1,2)=0D0
0186 P(1,3)=-Ebeam
0187 P(2,1)=0D0
0188 P(2,2)=0D0
0189 P(2,3)=Pbeam
0190
0191
0192 if(LQGPAR1(3).gt.0)then
0193 pxtot=sngl(P(1,1)+P(2,1))
0194 pytot=sngl(P(1,2)+P(2,2))
0195 pztot=sngl(P(1,3)+P(2,3))
0196 etot=sngl(dabs(P(1,3))+dabs(P(2,3)))
0197 ecmtot=sqrt(etot*etot-(pxtot*pxtot+pytot*pytot+pztot*pztot))
0198 if(pxtot.gt.0)then
0199 pxlow=pxtot-0.01*pxtot
0200 pxhi=pxtot+0.01*pxtot
0201 else
0202 pxlow=pxtot-1.
0203 pxhi=pxtot+1.
0204 endif
0205 if(pytot.gt.0)then
0206 pylow=pytot-0.01*pytot
0207 pyhi=pytot+0.01*pytot
0208 else
0209 pylow=pytot-1.
0210 pyhi=pytot+1.
0211 endif
0212 if(pztot.gt.0)then
0213 pzlow=pztot-0.01*pztot
0214 pzhi=pztot+0.01*pztot
0215 else
0216 pzlow=pztot-1.
0217 pzhi=pztot+1.
0218 endif
0219 if(etot.gt.0)then
0220 elow=etot-0.01*etot
0221 ehi=etot+0.01*etot
0222 else
0223 elow=etot-1.
0224 ehi=etot+1.
0225 endif
0226 if(ecmtot.gt.0)then
0227 ecmlow=ecmtot-0.01*ecmtot
0228 ecmhi=ecmtot+0.01*ecmtot
0229 else
0230 ecmlow=ecmtot-1.
0231 ecmhi=ecmtot+1.
0232 endif
0233 endif
0234
0235 if(LQGPAR1(5).eq.0)then
0236 Isub=401
0237 else
0238 Isub=LQGPAR1(5)
0239 endif
0240 sch=0.d0
0241 uch=0.d0
0242 call vzero(gproc,8)
0243 sigmax=LQGPAR3(12)
0244 if(echar.gt.0)then
0245 ibea=1
0246 else
0247 ibea=2
0248 endif
0249
0250
0251 CALL PYUPIN(ISUB,
0252 > CHBEA(ibea)//LQCHA(LQTYPE)//CHPROD(ibea,l_o),sigmax)
0253 MSEL=0
0254 MSUB(ISUB)=1
0255
0256 if(beampar(1).GT.0)then
0257 CALL PYINIT('USER','e+','p',0D0)
0258 else
0259 CALL PYINIT('USER','e-','p',0D0)
0260 endif
0261
0262
0263 if(LQGPAR1(3).ne.0)then
0264
0265 call hropen(69,'lqgenep','lqgenep.histo','N',1024,ierr)
0266 CALL hbook1(1000,'x gen',50,sngl(xmin),sngl(xmax),0.)
0267 CALL hbook1(1001,'x gen s-ch.',50,sngl(xmin),sngl(xmax),0.)
0268 CALL hbook1(1002,'x gen u-ch',50,sngl(xmin),sngl(xmax),0.)
0269 CALL hbook1(2000,'y gen',50,sngl(ymin),sngl(ymax),0.)
0270 CALL hbook1(2001,'y gen s-ch',50,sngl(ymin),sngl(ymax),0.)
0271 CALL hbook1(2002,'y gen u-ch',50,sngl(ymin),sngl(ymax),0.)
0272 CALL hbook1(3000,'Q2 gen',50,0.,6.,0.)
0273 call hbook1(5001,'sum px',100,pxlow,pxhi,0.)
0274 call hbook1(5002,'sum py',100,pylow,pyhi,0.)
0275 call hbook1(5003,'sum pz',100,pzlow,pzhi,0.)
0276 call hbook1(5004,'sum e',100,elow,ehi,0.)
0277 call hbook1(5000,'center of mass energy',
0278 > 100,ecmlow,ecmhi,0.)
0279 endif
0280 write(6,*)
0281 endif
0282
0283
0284 OPEN(asciiN, file=out_file)
0285
0286
0287 Integer HeaderTest = 0
0288 if(HeaderTest .LE. 1) then
0289 WRITE(29,*)' Pythia Output File'
0290 WRITE(29,*)'============================================'
0291 WRITE(29,*)'I, ievent, genevent, subprocess, nucleon,
0292 & targetparton, xtargparton, beamparton, xbeamparton,
0293 & thetabeamprtn, truey, trueQ2, truex, trueW2, trueNu,
0294 & leptonphi, s_hat, t_hat, u_hat, pt2_hat, Q2_hat, F2, F1, R,
0295 & sigma_rad, SigRadCor, EBrems, photonflux, nrTracks'
0296 WRITE(29,*)'============================================'
0297
0298
0299
0300
0301
0302 WRITE(29,30) 'I', 'KS', 'KF(ID)', 'ORIG', 'D1', 'D2', 'px'
0303 & , 'py', 'pz', 'E', 'm', 'vx', 'vy', 'vz'
0304 30 FORMAT(A5, A5, A7, A5, A5, A5, A10, A10, A10, A10, A10,
0305 & A10, A10, A10)
0306 WRITE(29,*)'============================================'
0307 endif
0308 HeaderTest = 999
0309
0310
0311
0312 if(flag.eq.1)Then
0313 CALL PYEVNT
0314
0315 CALL PYHEPC(1)
0316
0317
0318 if(genproc.eq.1)sch=sch+1
0319 if(genproc.eq.2)uch=uch+1
0320
0321
0322 gproc(genprtyp)=gproc(genprtyp)+1
0323
0324
0325 LQGKPAR(1)=X
0326 LQGKPAR(2)=Y
0327 LQGKPAR(3)=Q2
0328 LQGDDX(1)=(DXSEC(2)+DXSEC(3))*1.d-9
0329 LQGDDX(2)=DXSEC(3)*1.d-9
0330 LQGDDX(3)=DXSEC(2)*1.d-9
0331 LQGPID(1)=q_s
0332 LQGPID(2)=q_o
0333 if(genproc.eq.1)then
0334 LQGPID(3)=1
0335 elseif(genproc.eq.2)then
0336 LQGPID(3)=2
0337 endif
0338
0339
0340
0341
0342
0343 do 60 J=1,N
0344
0345 if ((K(J,1).EQ.11).and.
0346 > (K(J,2).EQ.15)) then ! find tau and get it's line number
0347 idt=J
0348 idtd=K(J,4) ! tau decay's to what line number
0349 endif
0350
0351 do 45 l=1,N
0352 if (l.EQ.idt) then
0353
0354 ptt=PYP(l,10)
0355 pxt=P(l,1)
0356 pyt=P(l,2)
0357 pht=PYP(l,16)
0358 ppt=PYP(l,14)
0359 endif
0360 45 enddo
0361
0362 if ((K(J,1).EQ.1).and.
0363 > (K(J,2).EQ.16)
0364 > .and.(K(J,3).EQ.idt)) then ! find tau neutrino and get it's line number
0365 idtnu=J
0366 endif
0367
0368 if ((J.EQ.idtd) ! line number is the decay of tau
0369 > .and.(K(J,2).NE.-12) ! this decay is not nu_ebar
0370 > .and.(K(J,2).NE.-14)) then ! this decay is not nu_mubar
0371
0372 ptid=0.d0
0373 pxid=0.d0
0374 pyid=0.d0
0375 pzid=0.d0
0376
0377 do 50 I=1,N
0378
0379 if ((K(I,1).LT.11) ! Anything that doesn't decay = Final Product
0380 > .and.(K(I,2).NE.12) ! Final Product NOT AN Electron neutrino
0381 > .and.(K(I,2).NE.14) ! Final Product NOT AN Muon neutrino
0382 > .and.(K(I,2).NE.16) ! Final Product NOT AN Tau neutrino
0383 > .and.(K(I,2).NE.-12) ! Final Product NOT AN Electron neutrino
0384 > .and.(K(I,2).NE.-14) ! Final Product NOT AN Muon neutrino
0385 > .and.(K(I,2).NE.-16)) then ! Final Product NOT A Tau neutrino
0386
0387 ptid=ptid+PYP(I,10)
0388 pxid=pxid+P(I,1)
0389 pyid=pyid+P(I,2)
0390 pzid=pzid+P(I,3)
0391 endif
0392 50 enddo
0393
0394 !!Delta Phi < 20 cut
0395
0396 if (pyid.GT.0.) then
0397 phimiss=acos(pxid/sqrt(pxid*pxid+pyid*pyid))*180/(3.14159)
0398 endif
0399
0400 if (pyid.LT.0.) then
0401 phimiss=-acos(pxid/sqrt(pxid*pxid+pyid*pyid))*180/(3.14159)
0402
0403 endif
0404
0405 if(abs(phimiss-pht).GT.180.) then
0406 dphimiss= 360-abs(phimiss-pht)
0407 endif
0408
0409 if(abs(phimiss-pht).GT.180.) then
0410 dphimiss= abs(phimiss-pht)
0411 endif
0412
0413
0414
0415 if (dphimiss.LT.20.) then
0416 write(8,10)LQGPAR1(4),pxid,pyid,pxt,
0417 > pyt,pht,ppt
0418 endif
0419
0420 10 FORMAT(I8,6(1PE14.6))
0421 endif
0422 60 enddo
0423
0424
0425 LQGPAR1(4)=LQGPAR1(4)+1
0426 if(LQGPAR1(4).LE.LQGPAR1(2)) CALL PYLIST(2)
0427
0428 if(mod(LQGPAR1(4),1000).eq.0)then
0429 write(6,1000) LQGPAR1(4)
0430 1000 format('>>>>>> ',I8,
0431 > ' events succesfully generated <<<<<<')
0432 endif
0433
0434
0435
0436
0437
0438
0439 WRITE(29,*) ' 0 1 1 95 2212 21 0.000000
0440 & 21 0.000000 0.000000 0.03587194790
0441 & 0.00016131040 0.00000022484 718.31915283203
0442 & 382.32009887695
0443 & 0.00773132309 0.00000000000 0.000000000
0444 & 0.000000000 0.000000000 0.000000000
0445 & 0.000000000 0.000000000 0.000000000
0446 & 0.000000000 0.000000000 0.000000000
0447 & 1.422321899 670'
0448 WRITE(29,*)'============================================'
0449 DO 100 L=1,N
0450
0451
0452 if(L == 3) then
0453 WRITE(29,*)'============================================'
0454 endif
0455
0456
0457
0458
0459
0460
0461 WRITE(29,40) L, K(L,1), K(L,2), K(L,3), K(L,4), K(L,5),
0462 & P(L,1), P(L,2)
0463 & , P(L,3), P(L,4), P(L,5), V(L,1), V(L,2), V(L,3)
0464 40 FORMAT(I5, I5, I7, I5, I5, I5, F10.5, F10.5, F10.5
0465 & , F10.5, F10.5, F15.5, F15.5, F15.5)
0466 KSprev = K(L,1)
0467 100 continue
0468 WRITE(29,*)'=============== Event finished ==============='
0469 CALL PYHIST
0470
0471 if(LQGPAR1(3).ne.0)then
0472
0473 CALL HF1(1000,sngl(x),1.)
0474 if(genproc.eq.1)CALL HF1(1001,sngl(x),1.)
0475 if(genproc.eq.2)CALL HF1(1002,sngl(x),1.)
0476 CALL HF1(2000,sngl(y),1.)
0477 if(genproc.eq.1)CALL HF1(2001,sngl(y),1.)
0478 if(genproc.eq.2)CALL HF1(2002,sngl(y),1.)
0479 CALL HF1(3000,log10(sngl(Q2)),1.)
0480
0481
0482 px_sum=0.
0483 py_sum=0.
0484 pz_sum=0.
0485 e_sum=0.
0486 cme=0.
0487 do 222 i=1,N
0488 if(K(I,1).le.10)then
0489 px_sum=px_sum+P(I,1)
0490 py_sum=py_sum+P(I,2)
0491 pz_sum=pz_sum+P(I,3)
0492 e_sum=e_sum+P(I,4)
0493 endif
0494 222 enddo
0495 cme=sqrt(e_sum**2-px_sum**2-py_sum**2-pz_sum**2)
0496 call hf1(5001,sngl(px_sum),1.)
0497 call hf1(5002,sngl(py_sum),1.)
0498 call hf1(5003,sngl(pz_sum),1.)
0499 call hf1(5004,sngl(e_sum),1.)
0500 call hf1(5000,sngl(cme),1.)
0501 endif
0502 endif
0503
0504
0505
0506 if(flag.eq.2)Then
0507 write(6,*)
0508
0509
0510 CALL PYSTAT(1)
0511
0512 write(6,*)
0513
0514 CALL LQGPRI2
0515
0516
0517 if(LQGPAR1(3).ne.0)then
0518 Call HCDIR('//lqgenep',' ')
0519 CALL HROUT(0,ICYCLE,' ')
0520 CALL HREND('lqgenep')
0521 endif
0522 endif
0523 END
0524
0525
0526
0527 SUBROUTINE PYUPEV(ISUB,SIGEV)
0528
0529
0530
0531
0532 IMPLICIT DOUBLE PRECISION(A-H, O-Z)
0533 IMPLICIT INTEGER(I-N)
0534 INTEGER PYK,PYCHGE,PYCOMP
0535
0536
0537 integer echar
0538 double precision ebeam,pbeam
0539 common /LQGbeam/ebeam,pbeam,echar
0540
0541
0542 double precision pdfsf(20)
0543 common /LQGpdfC/ pdfsf
0544
0545
0546 double precision xmax,xmin,ymax,ymin,zmax,zmin,Q2min
0547 common /LQGKinC/ xmax,xmin,ymax,ymin,zmax,zmin,Q2min
0548
0549
0550 double precision Mlq,G1,G2
0551 Integer LQtype,l_o,q_i,q_j
0552 common /LQGproc/ Mlq,G1,G2,LQtype,l_o,q_i,q_j
0553
0554
0555 double precision S,Srad,x,y,z,Q2
0556 common /LQGKinV/ S,Srad,x,y,z,Q2
0557
0558
0559 double precision DXSEC(3),pvalence
0560 integer q_o,q_s,genproc,genprtyp,sch,uch,gproc(8)
0561 common /LQGout/ DXSEC,pvalence,q_o,q_s,genproc,genprtyp,
0562 >sch,uch,gproc
0563
0564 CHARACTER CHAF*16
0565 COMMON /PYDAT4/CHAF(500,2)
0566 COMMON/PYDAT1/MSTU(200),PARU(200),MSTJ(200),PARJ(200)
0567 COMMON/PYUPPR/NUP,KUP(20,7),NFUP,IFUP(10,2),PUP(20,5),Q2UP(0:10)
0568 SAVE /PYDAT1/,/PYUPPR/
0569
0570 DIMENSION XPPs(-25:25),XPPu(-25:25),XPE(-25:25),TERM(20)
0571 DATA PI/3.141592653589793D0/
0572
0573 DATA CONV/1.D-9/
0574
0575
0576
0577 double precision BEAMPAR,LQGPAR3
0578 integer LQGPAR1,LQGPAR2
0579 CHARACTER*256 out_file
0580 COMMON/LQGDATA/BEAMPAR(3),LQGPAR1(10),LQGPAR2(10),LQGPAR3(20),
0581 > out_file
0582
0583
0584 Character*7 LQCHA(14)
0585 DATA LQCHA /'S_0L','S_0R','~S_0R','S_1L',
0586 > 'V_1/2L','V_1/2R','~V_1/2L',
0587 > 'V_0L','V_0R','~V_0R','V_1L',
0588 > 'S_1/2L','S_1/2R','~S_1/2L'/
0589
0590 DATA KLQ /39/
0591
0592 sigev=0.d0
0593 irej=0
0594
0595
0596 X=pyr(0)*(Xmax-Xmin)+Xmin
0597 Y=pyr(0)*(Ymax-Ymin)+Ymin
0598 Z=1
0599 Srad=S*z
0600 Q2=S*X*Y*Z
0601
0602
0603 call LQGDDXS
0604
0605 dxdydz=(Xmax-Xmin)*(Ymax-Ymin)*(Zmax-Zmin)
0606 Sigev=(DXSEC(2)+DXSEC(3))*conv*dxdydz
0607
0608
0609
0610 ECM2XZ=S*X*Z
0611 ECMXZ=sqrt(ECM2XZ)
0612 NUP=5
0613 KUP(1,1)=1
0614 KUP(1,2)=(echar)*-11
0615 KUP(1,3)=0
0616 KUP(1,4)=0
0617 KUP(1,5)=0
0618 KUP(1,6)=0
0619 KUP(1,7)=0
0620 PUP(1,1)=0.
0621 PUP(1,2)=0.
0622 PUP(1,4)=Z*sqrt(S)/2.d0
0623 PUP(1,3)=PUP(1,4)
0624 PUP(1,5)=0.
0625
0626 KUP(2,1)=1
0627 KUP(2,2)=q_s
0628 KUP(2,3)=0
0629 KUP(2,4)=0
0630 KUP(2,5)=0
0631 KUP(2,6)=0
0632 KUP(2,7)=0
0633 if(q_s.gt.0)then
0634 KUP(2,6)=5
0635 else
0636 KUP(2,7)=5
0637 endif
0638 PUP(2,1)=0.
0639 PUP(2,2)=0.
0640 PUP(2,4)=X*sqrt(S)/2.d0
0641 PUP(2,3)=-PUP(2,4)
0642 PUP(2,5)=0.
0643
0644 KUP(3,1)=2
0645 KUP(3,2)=KLQ
0646 CHAF(pycomp(KLQ),1)=LQCHA(LQTYPE)
0647 KUP(3,3)=0
0648 KUP(3,4)=0
0649 KUP(3,5)=0
0650 KUP(3,6)=0
0651 KUP(3,7)=0
0652 PUP(3,1)=PUP(2,1)+PUP(1,1)
0653 PUP(3,2)=PUP(2,2)+PUP(1,2)
0654 PUP(3,3)=PUP(2,3)+PUP(1,3)
0655 PUP(3,4)=sqrt(ECM2XZ+PUP(3,1)**2+PUP(3,2)**2+PUP(3,3)**2)
0656 PUP(3,5)=sqrt(ECM2XZ)
0657
0658
0659
0660 theta=acos(1.d0-2.d0*Y)
0661 PHI=2D0*PI*PYR(0)
0662 rtshat=ECMXZ
0663 KUP(4,1)=1
0664 KUP(4,2)=echar*-(11+2*(l_o-1))
0665 KUP(5,1)=1
0666 KUP(5,2)=q_o
0667 PUP(4,5)=PYMASS(KUP(4,2))
0668 PUP(5,5)=PYMASS(KUP(5,2))
0669 PUP44=0.5D0*(RTSHAT**2+PUP(4,5)**2-PUP(5,5)**2)/RTSHAT
0670 PUP54=RTSHAT-PUP44
0671 KUP(4,3)=3
0672 KUP(4,4)=0
0673 KUP(4,5)=0
0674 KUP(4,6)=0
0675 KUP(4,7)=0
0676 if(irej.eq.1.and.PUP44**2-PUP(4,5)**2.lt.0)then
0677 PMOD=1.d0
0678 else
0679 PMOD=sqrt(PUP44**2-PUP(4,5)**2)
0680 endif
0681 PUP(4,1)=PMOD*sin(theta)*cos(phi)
0682 PUP(4,2)=PMOD*sin(theta)*sin(phi)
0683 PUP43=PMOD*cos(theta)
0684 PUP44=PUP(3,5)/2.d0
0685
0686 KUP(5,3)=3
0687 KUP(5,4)=0
0688 KUP(5,5)=0
0689 KUP(5,6)=0
0690 KUP(5,7)=0
0691 if(q_o.gt.0)then
0692 KUP(5,4)=2
0693 else
0694 KUP(5,5)=2
0695 endif
0696 PUP(5,1)=0.
0697 PUP(5,2)=0.
0698 if(irej.eq.1.and.PUP54**2-PUP(5,5)**2.lt.0)then
0699 PMOD=1.d0
0700 else
0701 PMOD=sqrt(PUP54**2-PUP(5,5)**2)
0702 endif
0703 PUP(5,1)=-PUP(4,1)
0704 PUP(5,2)=-PUP(4,2)
0705 PUP53=-PUP43
0706
0707 beta=(z-x)/(z+x)
0708 gamma=0.5d0*(z+x)/sqrt(x*z)
0709 PUP(4,3)=GAMMA*(PUP43+BETA*PUP44)
0710 PUP(4,4)=GAMMA*(PUP44+BETA*PUP43)
0711 PUP(5,3)=GAMMA*(PUP53+BETA*PUP54)
0712 PUP(5,4)=GAMMA*(PUP54+BETA*PUP53)
0713
0714 NFUP=1
0715 IFUP(1,1)=4
0716 IFUP(1,2)=5
0717 Q2UP(0)=Q2
0718 Q2UP(1)=Q2
0719 RETURN
0720 END
0721
0722
0723 SUBROUTINE LQGDDXS
0724
0725
0726
0727
0728
0729 implicit none
0730
0731
0732 integer echar
0733 double precision ebeam,pbeam
0734 common /LQGbeam/ebeam,pbeam,echar
0735
0736
0737 double precision pdfsf(20)
0738 common /LQGpdfC/ pdfsf
0739
0740
0741 double precision xmax,xmin,ymax,ymin,zmax,zmin,Q2min
0742 common /LQGKinC/ xmax,xmin,ymax,ymin,zmax,zmin,Q2min
0743
0744
0745 double precision Mlq,G1,G2
0746 Integer LQtype,l_o,q_i,q_j
0747 common /LQGproc/ Mlq,G1,G2,LQtype,l_o,q_i,q_j
0748
0749
0750 double precision S,Srad,x,y,z,Q2
0751 common /LQGKinV/ S,Srad,x,y,z,Q2
0752
0753
0754 double precision DXSEC(3),pvalence
0755 integer q_o,q_s,genproc,genprtyp,sch,uch,gproc(8)
0756 common /LQGout/ DXSEC,pvalence,q_o,q_s,genproc,genprtyp,
0757 >sch,uch,gproc
0758
0759 double precision DSIGMADXDY(4)
0760 double precision rand(1),sfrac1,sfrac2,ufrac1,ufrac2
0761 double precision pvalences_u,pvalences_d,pvalenceu_u,pvalenceu_d
0762
0763
0764
0765
0766
0767
0768
0769
0770
0771
0772
0773
0774
0775
0776
0777
0778
0779
0780
0781
0782
0783
0784
0785
0786
0787
0788
0789
0790
0791
0792
0793
0794
0795
0796
0797
0798
0799
0800
0801
0802
0803
0804
0805
0806
0807
0808
0809
0810
0811
0812
0813
0814 double precision pyr
0815 double precision C_R_P,C_L_P,C_R_E,C_L_E
0816 & ,C_R_U,C_L_U,C_R_D,C_L_D
0817 & ,B_RR_U,B_RL_U
0818 & ,B_LR_U,B_LL_U
0819 & ,B_RR_D,B_RL_D
0820 & ,B_LR_D,B_LL_D
0821 double precision CCCf2u,DDDf2u,EEEf2u,FFFf2u,
0822 & CCCf2d,DDDf2d,EEEf2d,FFFf2d
0823 double precision CCCf0u,DDDf0u,EEEf0u,FFFf0u,
0824 & CCCf0d,DDDf0d,EEEf0d,FFFf0d
0825 double precision CCCv2u,DDDv2u,EEEv2u,FFFv2u,
0826 & CCCv2d,DDDv2d,EEEv2d,FFFv2d
0827 double precision CCCv0u,DDDv0u,EEEv0u,FFFv0u,
0828 & CCCv0d,DDDv0d,EEEv0d,FFFv0d
0829
0830 double precision weakmix,Mz2,Gz2,Gf,alpha,A,P,pi
0831 parameter (weakmix=0.2315, Mz2=(91.187)**2,Gz2=(2.490)**2)
0832 parameter (Gf=0.0000116639,alpha=1.0/137.036, A=27.5, P=820.0)
0833 parameter (pi=3.141592653589793d0)
0834
0835 double precision UPV,DNV,USEA,DSEA,STR,CHM,BOT,TOP,GL
0836 double precision UPQs(3),DNQs(3),UPQBs(3),DNQBs(3)
0837 double precision UPQu(3),DNQu(3),UPQBu(3),DNQBu(3)
0838 double precision scales,scaleu,LambdaL2,LambdaR2,Ms2
0839 & ,aaa,bbb,ggg
0840 & ,LambdaL2_1,LambdaL2_2
0841 & ,LambdaR2_1,LambdaR2_2
0842 & ,GAM
0843
0844
0845 integer LRindex(14)
0846 data LRindex/0,1,1,0,0,1,0,0,1,1,0,0,1,0/
0847 integer SVindex(14)
0848 data SVindex/1,1,1,1,2,2,2,3,3,3,3,4,4,4/
0849
0850
0851 if(y.eq.1)y=1.d0-1.d-13
0852
0853 DSIGMADXDY(1)=0.d0
0854 DSIGMADXDY(2)=0.d0
0855 DSIGMADXDY(3)=0.d0
0856 DSIGMADXDY(4)=0.d0
0857 q_o=0
0858 pvalence=0.
0859 pvalences_u=0.
0860 pvalenceu_u=0.
0861 pvalences_d=0.
0862 pvalenceu_d=0.
0863
0864 C_R_P = weakmix
0865 C_L_P = -0.5+weakmix
0866 C_R_U = -2.0*weakmix/3.0
0867 C_L_U = 0.5-2.0*weakmix/3.0
0868 C_R_D = weakmix/3.0
0869 C_L_D = -0.5+weakmix/3.0
0870 C_R_E = weakmix
0871 C_L_E = -0.5+weakmix
0872
0873 Ms2=(MLQ)**2
0874 Q2=Srad*X*Y
0875 if(Q2.lt.Q2min)goto 999
0876
0877 scaleu=sqrt(Srad*X*(1.-Y))
0878 CALL STRUCTM(X,scaleu,UPV
0879 & ,DNV,USEA,DSEA,STR,CHM,BOT,TOP,GL)
0880 if(UPV+USEA.gt.0)
0881 & pvalenceu_u=UPV/(UPV+USEA)
0882 if(DNV+DSEA.gt.0)
0883 & pvalenceu_d=DNV/(DNV+DSEA)
0884 if(echar.eq.1)then
0885
0886 UPQu(1)=UPV+USEA
0887 UPQBu(1)=USEA
0888 UPQu(2)=CHM
0889 UPQBu(2)=CHM
0890 UPQu(3)=TOP
0891 UPQBu(3)=TOP
0892 DNQu(1)=DNV+DSEA
0893 DNQBu(1)=DSEA
0894 DNQu(2)=STR
0895 DNQBu(2)=STR
0896 DNQu(3)=BOT
0897 DNQBu(3)=BOT
0898 elseif(echar.eq.-1)then
0899
0900 UPQu(1)=USEA
0901 UPQBu(1)=UPV+USEA
0902 UPQu(2)=CHM
0903 UPQBu(2)=CHM
0904 UPQu(3)=TOP
0905 UPQBu(3)=TOP
0906 DNQu(1)=DSEA
0907 DNQBu(1)=DNV+DSEA
0908 DNQu(2)=STR
0909 DNQBu(2)=STR
0910 DNQu(3)=BOT
0911 DNQBu(3)=BOT
0912 endif
0913
0914 scales=sqrt(Srad*X)
0915 CALL STRUCTM(X,scales,UPV
0916 & ,DNV,USEA,DSEA,STR,CHM,BOT,TOP,GL)
0917 if(UPV+USEA.gt.0)
0918 & pvalences_u=UPV/(UPV+USEA)
0919 if(DNV+DSEA.gt.0)
0920 & pvalences_d=DNV/(DNV+DSEA)
0921 if(echar.eq.1)then
0922
0923 UPQs(1)=UPV+USEA
0924 UPQBs(1)=USEA
0925 UPQs(2)=CHM
0926 UPQBs(2)=CHM
0927 UPQs(3)=TOP
0928 UPQBs(3)=TOP
0929 DNQs(1)=DNV+DSEA
0930 DNQBs(1)=DSEA
0931 DNQs(2)=STR
0932 DNQBs(2)=STR
0933 DNQs(3)=BOT
0934 DNQBs(3)=BOT
0935 elseif(echar.eq.-1)then
0936
0937 UPQs(1)=USEA
0938 UPQBs(1)=UPV+USEA
0939 UPQs(2)=CHM
0940 UPQBs(2)=CHM
0941 UPQs(3)=TOP
0942 UPQBs(3)=TOP
0943 DNQs(1)=DSEA
0944 DNQBs(1)=DNV+DSEA
0945 DNQs(2)=STR
0946 DNQBs(2)=STR
0947 DNQs(3)=BOT
0948 DNQBs(3)=BOT
0949 endif
0950
0951 aaa = Q2*(Q2+Mz2)/((Q2+Mz2)**2+Mz2*Gz2)
0952 bbb = sqrt(2.0)*Gf*Mz2/(pi*alpha)
0953
0954 B_RR_U = -2./3. + aaa*bbb*C_R_P*C_R_U
0955 B_RL_U = -2./3. + aaa*bbb*C_R_P*C_L_U
0956 B_LR_U = -2./3. + aaa*bbb*C_L_P*C_R_U
0957 B_LL_U = -2./3. + aaa*bbb*C_L_P*C_L_U
0958 B_RR_D = 1./3. + aaa*bbb*C_R_P*C_R_D
0959 B_RL_D = 1./3. + aaa*bbb*C_R_P*C_L_D
0960 B_LR_D = 1./3. + aaa*bbb*C_L_P*C_R_D
0961 B_LL_D = 1./3. + aaa*bbb*C_L_P*C_L_D
0962
0963 IF (SVindex(LQTYPE).EQ.1) THEN
0964
0965 CCCf2u=Q2*(1.-Y)**2*(UPV+USEA)
0966 & /(4.0*pi*alpha)
0967 DDDf2u=Q2*USEA/(4.0*pi*alpha)
0968 EEEf2u=(Q2-X*Srad)**2*Y**2*(UPQu(q_j))
0969 & /(64.0*(pi*alpha)**2)
0970 FFFf2u=Q2**2*UPQBs(q_i)/(64.0*(pi*alpha)**2)
0971
0972 CCCf2d=Q2*(1.-Y)**2*(DNV+DSEA)
0973 & /(4.0*pi*alpha)
0974 DDDf2d=Q2*DSEA/(4.0*pi*alpha)
0975 EEEf2d=(Q2-X*Srad)**2*Y**2*(DNQu(q_j))
0976 & /(64.0*(pi*alpha)**2)
0977 FFFf2d=Q2**2*DNQBs(q_i)/(64.0*(pi*alpha)**2)
0978
0979 ELSE IF (SVindex(LQTYPE).EQ.4) THEN
0980
0981 CCCf0u=Q2*(1.-Y)**2*USEA/(4.0*pi*alpha)
0982 DDDf0u=Q2*(UPV+USEA)/(4.0*pi*alpha)
0983 EEEf0u=(Q2-X*Srad)**2*Y**2*UPQBu(q_j)
0984 & /(64.0*(pi*alpha)**2)
0985 FFFf0u=Q2**2*(UPQs(q_i))/(64.0*(pi*alpha)**2)
0986
0987 CCCf0d=Q2*(1.-Y)**2*DSEA/(4.0*pi*alpha)
0988 DDDf0d=Q2*(DNV+DSEA)/(4.0*pi*alpha)
0989 EEEf0d=(Q2-X*Srad)**2*Y**2*DNQBu(q_j)
0990 & /(64.0*(pi*alpha)**2)
0991 FFFf0d=Q2**2*(DNQs(q_i))/(64.0*(pi*alpha)**2)
0992
0993 ELSE IF (SVindex(LQTYPE).EQ.2) THEN
0994
0995 CCCv2u=Q2*(1.-Y)**2*USEA/(2.0*pi*alpha)
0996 DDDv2u=Q2*(UPV+USEA)/(2.0*pi*alpha)
0997 EEEv2u=Q2**2*UPQBs(q_i)/(16.0*(pi*alpha)**2)
0998 FFFv2u=(Q2-X*Srad)**2*Y**2*(UPQu(q_j))
0999 & /(16.0*(pi*alpha)**2*(1.0-Y)**2)
1000
1001 CCCv2d=Q2*(1.-Y)**2*DSEA/(2.0*pi*alpha)
1002 DDDv2d=Q2*(DNV+DSEA)/(2.0*pi*alpha)
1003 EEEv2d=Q2**2*DNQBs(q_i)/(16.0*(pi*alpha)**2)
1004 FFFv2d=(Q2-X*Srad)**2*Y**2*(DNQu(q_j))
1005 & /(16.0*(pi*alpha)**2*(1.0-Y)**2)
1006
1007 ELSE IF (SVindex(LQTYPE).EQ.3) THEN
1008
1009 CCCv0u=Q2*(1.-Y)**2*(UPV+USEA)
1010 & /(2.0*pi*alpha)
1011 DDDv0u=Q2*USEA/(2.0*pi*alpha)
1012 EEEv0u=Q2**2*(UPQs(q_i))/(16.0*(pi*alpha)**2)
1013 FFFv0u=(Q2-X*Srad)**2*Y**2*UPQBu(q_j)
1014 & /(16.0*(pi*alpha)**2*(1.0-Y)**2)
1015
1016 CCCv0d=Q2*(1.-Y)**2*(DNV+DSEA)
1017 & /(2.0*pi*alpha)
1018 DDDv0d=Q2*DSEA/(2.0*pi*alpha)
1019 EEEv0d=Q2**2*(DNQs(q_i))/(16.0*(pi*alpha)**2)
1020 FFFv0d=(Q2-X*Srad)**2*Y**2*DNQBu(q_j)
1021 & /(16.0*(pi*alpha)**2*(1.0-Y)**2)
1022 ENDIF
1023
1024 DSIGMADXDY(1)=(B_RL_U**2+B_LR_U**2+
1025 & (B_RR_U**2+B_LL_U**2)*(1.-Y)**2)
1026 & *(UPV+USEA+CHM+TOP)+
1027 & (B_RL_D**2+B_LR_D**2+
1028 & (B_RR_D**2+B_LL_D**2)*(1.-Y)**2)
1029 & *(DNV+DSEA+STR+BOT)+
1030 & (B_RR_U**2+B_LL_U**2+
1031 & (B_RL_U**2+B_LR_U**2)*(1.-Y)**2)
1032 & *(USEA+CHM+TOP) +
1033 & (B_RR_D**2+B_LL_D**2+
1034 & (B_RL_D**2+B_LR_D**2)*(1.-Y)**2)
1035 & *(DSEA+STR+BOT)
1036
1037 IF (LRindex(LQTYPE).EQ.0) THEN
1038 LambdaL2_1=G1
1039 LambdaL2_2=G2
1040 if(G2.ne.0)then
1041 LambdaL2=LambdaL2_1*LambdaL2_2
1042 else
1043 LambdaL2=G1*G1
1044 endif
1045 LambdaR2_1=0.0
1046 LambdaR2_2=0.0
1047 LambdaR2=LambdaR2_1*LambdaR2_2
1048 ELSE IF (LRindex(LQTYPE).EQ.1) THEN
1049 LambdaR2_1=G1
1050 LambdaR2_2=G2
1051 if(G2.ne.0)then
1052 LambdaR2=LambdaR2_1*LambdaR2_2
1053 else
1054 LambdaR2=G1*G1
1055 endif
1056 LambdaL2_1=0.0
1057 LambdaL2_2=0.0
1058 LambdaL2=LambdaL2_1*LambdaL2_2
1059 ENDIF
1060
1061 IF (LQTYPE.EQ.1.or.LQTYPE.EQ.2) THEN
1062 GAM=MLQ*(sqrt(2.0)*LambdaL2_1+LambdaR2_1)**2
1063 if(l_o.gt.1)GAM=GAM+
1064 & MLQ*(sqrt(2.0)*LambdaL2_2+LambdaR2_2)**2
1065 AAA=(X*Srad-Ms2)**2+
1066 & Ms2*GAM**2/((16.0*pi)**2)
1067 BBB=LambdaR2*B_RR_U+LambdaL2*B_LL_U
1068 DSIGMADXDY(2)=BBB*CCCf2u/(Q2-X*Srad-Ms2)+
1069 & BBB*DDDf2u*(X*Srad-Ms2)/AAA
1070 GGG=(LambdaR2+LambdaL2)**2
1071 if(q_i.ne.3)then
1072 DSIGMADXDY(3)=EEEf2u*GGG/(Q2-X*Srad-Ms2)**2
1073 else
1074
1075 DSIGMADXDY(3)=0.d0
1076 endif
1077 if(q_j.ne.3)then
1078 DSIGMADXDY(4)=FFFf2u*GGG/AAA
1079 else
1080
1081 DSIGMADXDY(4)=0.d0
1082 endif
1083
1084 if(DSIGMADXDY(3)+DSIGMADXDY(4).eq.0)then
1085 goto 999
1086 endif
1087 sfrac1=DSIGMADXDY(4)/(DSIGMADXDY(4)+DSIGMADXDY(3))
1088
1089 rand(1)=pyr(0)
1090 if(rand(1).lt.sfrac1)then
1091
1092 genproc=1
1093 genprtyp=3
1094 if(q_i.eq.1)then
1095 if(echar.lt.0)pvalence=pvalences_u
1096 endif
1097 if(q_j.eq.1)q_o=-2
1098 if(q_j.eq.2)q_o=-4
1099 if(q_j.eq.3)q_o=-6
1100 if(q_i.eq.1)q_s=-2
1101 if(q_i.eq.2)q_s=-4
1102 if(q_i.eq.3)q_s=-6
1103 else
1104
1105 genproc=2
1106 genprtyp=5
1107 if(q_j.eq.1)then
1108 if(echar.gt.0)pvalence=pvalenceu_u
1109 endif
1110 if(q_i.eq.1)q_o=2
1111 if(q_i.eq.2)q_o=4
1112 if(q_i.eq.3)q_o=6
1113 if(q_j.eq.1)q_s=2
1114 if(q_j.eq.2)q_s=4
1115 if(q_j.eq.3)q_s=6
1116 endif
1117 ELSE IF (LQTYPE.EQ.3) THEN
1118 GAM=MLQ*(LambdaL2_1+LambdaR2_1)**2
1119 if(l_o.gt.1)GAM=GAM+
1120 & MLQ*(LambdaL2_2+LambdaR2_2)**2
1121 AAA=(X*Srad-Ms2)**2+
1122 & Ms2*GAM**2/((16.0*pi)**2)
1123 BBB=LambdaR2*B_RR_D+LambdaL2*B_LL_D
1124 DSIGMADXDY(2)=BBB*CCCf2d/(Q2-X*Srad-Ms2)+
1125 & BBB*DDDf2d*(X*Srad-Ms2)/AAA
1126 GGG=(LambdaR2+LambdaL2)**2
1127 DSIGMADXDY(3)=EEEf2d*GGG/(Q2-X*Srad-Ms2)**2
1128 DSIGMADXDY(4)=FFFf2d*GGG/AAA
1129
1130 if(DSIGMADXDY(3)+DSIGMADXDY(4).eq.0)then
1131 goto 999
1132 endif
1133 sfrac1=DSIGMADXDY(4)/(DSIGMADXDY(4)+DSIGMADXDY(3))
1134
1135 rand(1)=pyr(0)
1136 if(rand(1).lt.sfrac1)then
1137
1138 genproc=1
1139 genprtyp=4
1140 if(q_i.eq.1)then
1141 if(echar.lt.0)pvalence=pvalences_d
1142 endif
1143 if(q_j.eq.1)q_o=-1
1144 if(q_j.eq.2)q_o=-3
1145 if(q_j.eq.3)q_o=-5
1146 if(q_i.eq.1)q_s=-1
1147 if(q_i.eq.2)q_s=-3
1148 if(q_i.eq.3)q_s=-5
1149 else
1150
1151 genproc=2
1152 genprtyp=6
1153 if(q_j.eq.1)then
1154 if(echar.gt.0)pvalence=pvalenceu_d
1155 endif
1156 if(q_i.eq.1)q_o=1
1157 if(q_i.eq.2)q_o=3
1158 if(q_i.eq.3)q_o=5
1159 if(q_j.eq.1)q_s=1
1160 if(q_j.eq.2)q_s=3
1161 if(q_j.eq.3)q_s=5
1162 endif
1163 ELSE IF (LQTYPE.EQ.4) THEN
1164 GAM=MLQ*(sqrt(2.0)*LambdaL2_1+LambdaR2_1)**2
1165 if(l_o.gt.1)GAM=GAM+
1166 & MLQ*(sqrt(2.0)*LambdaL2_2+LambdaR2_2)**2
1167 AAA=(X*Srad-Ms2)**2+
1168 & Ms2*GAM**2/((16.0*pi)**2)
1169 BBB=LambdaR2*B_RR_D+2.0*LambdaL2*B_LL_D
1170 DSIGMADXDY(2)=BBB*CCCf2d/(Q2-X*Srad-Ms2)+
1171 & BBB*DDDf2d*(X*Srad-Ms2)/AAA
1172 GGG=(LambdaR2+2.0*LambdaL2)**2
1173 DSIGMADXDY(3)=EEEf2d*GGG/(Q2-X*Srad-Ms2)**2
1174 DSIGMADXDY(4)=FFFf2d*GGG/AAA
1175 sfrac1=DSIGMADXDY(4)
1176 ufrac1=DSIGMADXDY(3)
1177 AAA=(X*Srad-Ms2)**2+
1178 & Ms2*GAM**2/((16.0*pi)**2)
1179 BBB=LambdaR2*B_RR_U+LambdaL2*B_LL_U
1180 DSIGMADXDY(2)=DSIGMADXDY(2)+BBB*CCCf2u/(Q2-X*Srad-Ms2)+
1181 & BBB*DDDf2u*(X*Srad-Ms2)/AAA
1182 GGG=(LambdaR2+LambdaL2)**2
1183
1184 if(q_i.ne.3)then
1185 DSIGMADXDY(3)=DSIGMADXDY(3)+EEEf2u*GGG/(Q2-X*Srad-Ms2)**2
1186 endif
1187
1188 if(q_j.ne.3)then
1189 DSIGMADXDY(4)=DSIGMADXDY(4)+FFFf2u*GGG/AAA
1190 endif
1191
1192 if(DSIGMADXDY(3)+DSIGMADXDY(4).eq.0)then
1193 goto 999
1194 endif
1195 sfrac2=DSIGMADXDY(4)-sfrac1
1196 ufrac2=DSIGMADXDY(3)-ufrac1
1197 sfrac1=sfrac1/(DSIGMADXDY(3)+DSIGMADXDY(4))
1198 sfrac2=sfrac2/(DSIGMADXDY(3)+DSIGMADXDY(4))
1199 ufrac1=ufrac1/(DSIGMADXDY(3)+DSIGMADXDY(4))
1200 ufrac2=ufrac2/(DSIGMADXDY(3)+DSIGMADXDY(4))
1201
1202 rand(1)=pyr(0)
1203 if(rand(1).lt.sfrac1)then
1204
1205 genproc=1
1206 genprtyp=4
1207 if(q_i.eq.1)then
1208 if(echar.lt.0)pvalence=pvalences_d
1209 endif
1210 if(q_j.eq.1)q_o=-1
1211 if(q_j.eq.2)q_o=-3
1212 if(q_j.eq.3)q_o=-5
1213 if(q_i.eq.1)q_s=-1
1214 if(q_i.eq.2)q_s=-3
1215 if(q_i.eq.3)q_s=-5
1216 elseif(rand(1).lt.(sfrac1+sfrac2))then
1217
1218 genproc=1
1219 genprtyp=3
1220 if(q_i.eq.1)then
1221 if(echar.lt.0)pvalence=pvalences_u
1222 endif
1223 if(q_j.eq.1)q_o=-2
1224 if(q_j.eq.2)q_o=-4
1225 if(q_j.eq.3)q_o=-6
1226 if(q_i.eq.1)q_s=-2
1227 if(q_i.eq.2)q_s=-4
1228 if(q_i.eq.3)q_s=-6
1229 elseif(rand(1).lt.(sfrac1+sfrac2+ufrac1))then
1230
1231 genproc=2
1232 genprtyp=6
1233 if(q_j.eq.1)then
1234 if(echar.gt.0)pvalence=pvalenceu_d
1235 endif
1236 if(q_i.eq.1)q_o=1
1237 if(q_i.eq.2)q_o=3
1238 if(q_i.eq.3)q_o=5
1239 if(q_j.eq.1)q_s=1
1240 if(q_j.eq.2)q_s=3
1241 if(q_j.eq.3)q_s=5
1242 else
1243
1244 genproc=2
1245 genprtyp=5
1246 if(q_j.eq.1)then
1247 if(echar.gt.0)pvalence=pvalenceu_u
1248 endif
1249 if(q_i.eq.1)q_o=2
1250 if(q_i.eq.2)q_o=4
1251 if(q_i.eq.3)q_o=6
1252 if(q_j.eq.1)q_s=2
1253 if(q_j.eq.2)q_s=4
1254 if(q_j.eq.3)q_s=6
1255 endif
1256 ELSE IF (LQTYPE.EQ.5) THEN
1257 AAA=(X*Srad-Ms2)
1258 GAM=MLQ*(LambdaL2_1+LambdaR2_1)**2
1259 if(l_o.gt.1)GAM=GAM+
1260 & MLQ*(LambdaL2_2+LambdaR2_2)**2
1261 BBB=LambdaR2*B_RL_D+LambdaL2*B_LR_D
1262 GGG=(LambdaL2+LambdaR2)/(24.0*pi)
1263 DSIGMADXDY(2)=BBB*CCCv2d*AAA/(AAA**2+
1264 & (Ms2*GGG)**2)+
1265 & BBB*DDDv2d/(Q2-X*Srad-Ms2)
1266 DSIGMADXDY(4)=EEEv2d*
1267 & ((LambdaL2**2+LambdaR2**2)*(1.-Y)**2+
1268 & 2.0*LambdaL2*LambdaR2*Y**2)/
1269 & (AAA**2+Ms2*GAM**2/(24.0*pi)**2)
1270 DSIGMADXDY(3)=FFFv2d*(LambdaL2**2+LambdaR2**2+
1271 & 2.0*Y**2*LambdaL2*LambdaR2)/
1272 & (Q2-X*Srad-Ms2)**2
1273
1274 if(DSIGMADXDY(3)+DSIGMADXDY(4).eq.0)then
1275 goto 999
1276 endif
1277 sfrac1=DSIGMADXDY(4)/(DSIGMADXDY(3)+DSIGMADXDY(4))
1278
1279 rand(1)=pyr(0)
1280 if(rand(1).lt.sfrac1)then
1281
1282 genproc=1
1283 genprtyp=4
1284 if(q_i.eq.1)then
1285 if(echar.lt.0)pvalence=pvalences_d
1286 endif
1287 if(q_j.eq.1)q_o=-1
1288 if(q_j.eq.2)q_o=-3
1289 if(q_j.eq.3)q_o=-5
1290 if(q_i.eq.1)q_s=-1
1291 if(q_i.eq.2)q_s=-3
1292 if(q_i.eq.3)q_s=-5
1293 else
1294
1295 genproc=2
1296 genprtyp=6
1297 if(q_j.eq.1)then
1298 if(echar.gt.0)pvalence=pvalenceu_d
1299 endif
1300 if(q_i.eq.1)q_o=1
1301 if(q_i.eq.2)q_o=3
1302 if(q_i.eq.3)q_o=5
1303 if(q_j.eq.1)q_s=1
1304 if(q_j.eq.2)q_s=3
1305 if(q_j.eq.3)q_s=5
1306 endif
1307 ELSE IF (LQTYPE.EQ.6) THEN
1308 AAA=(X*Srad-Ms2)
1309 GAM=MLQ*(LambdaL2_1+LambdaR2_1)**2
1310 if(l_o.gt.1)GAM=GAM+
1311 & MLQ*(LambdaL2_2+LambdaR2_2)**2
1312 BBB=LambdaR2*B_RL_D+LambdaL2*B_LR_D
1313 GGG=(LambdaL2+LambdaR2)/(24.0*pi)
1314 DSIGMADXDY(2)=BBB*CCCv2d*AAA/(AAA**2+
1315 & (Ms2*GGG)**2)+
1316 & BBB*DDDv2d/(Q2-X*Srad-Ms2)
1317 DSIGMADXDY(4)=EEEv2d*
1318 & ((LambdaL2**2+LambdaR2**2)*(1.-Y)**2+
1319 & 2.0*LambdaL2*LambdaR2*Y**2)/
1320 & (AAA**2+Ms2*GAM**2/(24.0*pi)**2)
1321 DSIGMADXDY(3)=FFFv2d*(LambdaL2**2+LambdaR2**2+
1322 & 2.0*Y**2*LambdaL2*LambdaR2)/
1323 & (Q2-X*Srad-Ms2)**2
1324
1325 sfrac1=DSIGMADXDY(4)
1326 ufrac1=DSIGMADXDY(3)
1327 BBB=LambdaR2*B_RL_U+LambdaL2*B_LR_U
1328 GGG=(LambdaL2+LambdaR2)/(24.0*pi)
1329 DSIGMADXDY(2)=DSIGMADXDY(2)+BBB*CCCv2u*AAA/(AAA**2+
1330 & (Ms2*GGG)**2)+
1331 & BBB*DDDv2u/(Q2-X*Srad-Ms2)
1332
1333 if(q_j.ne.3)then
1334 DSIGMADXDY(4)=DSIGMADXDY(4)+EEEv2u*
1335 & ((LambdaL2**2+LambdaR2**2)*(1.-Y)**2+
1336 & 2.0*LambdaL2*LambdaR2*Y**2)/
1337 & (AAA**2+Ms2*GAM**2/(24.0*pi)**2)
1338 endif
1339
1340 if(q_i.ne.3)then
1341 DSIGMADXDY(3)=DSIGMADXDY(3)+FFFv2u*(LambdaL2**2+LambdaR2**2+
1342 & 2.0*Y**2*LambdaL2*LambdaR2)/
1343 & (Q2-X*Srad-Ms2)**2
1344 endif
1345
1346 if(DSIGMADXDY(3)+DSIGMADXDY(4).eq.0)then
1347 goto 999
1348 endif
1349 sfrac2=DSIGMADXDY(4)-sfrac1
1350 ufrac2=DSIGMADXDY(3)-ufrac1
1351 sfrac1=sfrac1/(DSIGMADXDY(3)+DSIGMADXDY(4))
1352 sfrac2=sfrac2/(DSIGMADXDY(3)+DSIGMADXDY(4))
1353 ufrac1=ufrac1/(DSIGMADXDY(3)+DSIGMADXDY(4))
1354 ufrac2=ufrac2/(DSIGMADXDY(3)+DSIGMADXDY(4))
1355
1356 rand(1)=pyr(0)
1357 if(rand(1).lt.sfrac1)then
1358
1359 genproc=1
1360 genprtyp=4
1361 if(q_i.eq.1)then
1362 if(echar.lt.0)pvalence=pvalences_d
1363 endif
1364 if(q_j.eq.1)q_o=-1
1365 if(q_j.eq.2)q_o=-3
1366 if(q_j.eq.3)q_o=-5
1367 if(q_i.eq.1)q_s=-1
1368 if(q_i.eq.2)q_s=-3
1369 if(q_i.eq.3)q_s=-5
1370 elseif(rand(1).lt.(sfrac1+sfrac2))then
1371
1372 genproc=1
1373 genprtyp=3
1374 if(q_i.eq.1)then
1375 if(echar.lt.0)pvalence=pvalences_u
1376 endif
1377 if(q_j.eq.1)q_o=-2
1378 if(q_j.eq.2)q_o=-4
1379 if(q_j.eq.3)q_o=-6
1380 if(q_i.eq.1)q_s=-2
1381 if(q_i.eq.2)q_s=-4
1382 if(q_i.eq.3)q_s=-6
1383 elseif(rand(1).lt.(sfrac1+sfrac2+ufrac1))then
1384
1385 genproc=2
1386 genprtyp=6
1387 if(q_j.eq.1)then
1388 if(echar.gt.0)pvalence=pvalenceu_d
1389 endif
1390 if(q_i.eq.1)q_o=1
1391 if(q_i.eq.2)q_o=3
1392 if(q_i.eq.3)q_o=5
1393 if(q_j.eq.1)q_s=1
1394 if(q_j.eq.2)q_s=3
1395 if(q_j.eq.3)q_s=5
1396 else
1397
1398 genproc=2
1399 genprtyp=5
1400 if(q_j.eq.1)then
1401 if(echar.gt.0)pvalence=pvalenceu_u
1402 endif
1403 if(q_i.eq.1)q_o=2
1404 if(q_i.eq.2)q_o=4
1405 if(q_i.eq.3)q_o=6
1406 if(q_j.eq.1)q_s=2
1407 if(q_j.eq.2)q_s=4
1408 if(q_j.eq.3)q_s=6
1409 endif
1410 ELSE IF (LQTYPE.EQ.7) THEN
1411 AAA=(X*Srad-Ms2)
1412 GAM=MLQ*(LambdaL2_1+LambdaR2_1)**2
1413 if(l_o.gt.1)GAM=GAM+
1414 & MLQ*(LambdaL2_2+LambdaR2_2)**2
1415 BBB=LambdaR2*B_RL_U+LambdaL2*B_LR_U
1416 GGG=(LambdaL2+LambdaR2)/(24.0*pi)
1417 DSIGMADXDY(2)=BBB*CCCv2u*AAA/(AAA**2+
1418 & (Ms2*GGG)**2)+
1419 & BBB*DDDv2u/(Q2-X*Srad-Ms2)
1420
1421 if(q_j.ne.3)then
1422 DSIGMADXDY(4)=EEEv2u*
1423 & ((LambdaL2**2+LambdaR2**2)*(1.-Y)**2+
1424 & 2.0*LambdaL2*LambdaR2*Y**2)/
1425 & (AAA**2+Ms2*GAM**2/(24.0*pi)**2)
1426 else
1427 DSIGMADXDY(4)=0.d0
1428 endif
1429
1430 if(q_i.ne.3)then
1431 DSIGMADXDY(3)=FFFv2u*(LambdaL2**2+LambdaR2**2+
1432 & 2.0*Y**2*LambdaL2*LambdaR2)/
1433 & (Q2-X*Srad-Ms2)**2
1434 else
1435 DSIGMADXDY(3)=0.d0
1436 endif
1437
1438 if(DSIGMADXDY(3)+DSIGMADXDY(4).eq.0)then
1439 goto 999
1440 endif
1441 sfrac1=DSIGMADXDY(4)/(DSIGMADXDY(3)+DSIGMADXDY(4))
1442
1443 rand(1)=pyr(0)
1444 if(rand(1).lt.sfrac1)then
1445
1446 genproc=1
1447 genprtyp=3
1448 if(q_i.eq.1)then
1449 if(echar.lt.0)pvalence=pvalences_u
1450 endif
1451 if(q_j.eq.1)q_o=-2
1452 if(q_j.eq.2)q_o=-4
1453 if(q_j.eq.3)q_o=-6
1454 if(q_i.eq.1)q_s=-2
1455 if(q_i.eq.2)q_s=-4
1456 if(q_i.eq.3)q_s=-6
1457 else
1458
1459 genproc=2
1460 genprtyp=5
1461 if(q_j.eq.1)then
1462 if(echar.gt.0)pvalence=pvalenceu_u
1463 endif
1464 if(q_i.eq.1)q_o=2
1465 if(q_i.eq.2)q_o=4
1466 if(q_i.eq.3)q_o=6
1467 if(q_j.eq.1)q_s=2
1468 if(q_j.eq.2)q_s=4
1469 if(q_j.eq.3)q_s=6
1470 endif
1471 ELSE IF (LQTYPE.EQ.8.or.LQTYPE.EQ.9) THEN
1472 AAA=(X*Srad-Ms2)
1473 GAM=MLQ*(sqrt(2.0)*LambdaL2_1+LambdaR2_1)**2
1474 if(l_o.gt.1)GAM=GAM+
1475 & MLQ*(sqrt(2.0)*LambdaL2_2+LambdaR2_2)**2
1476 BBB=-LambdaR2*B_RR_D-LambdaL2*B_LL_D
1477 GGG=(2.0*LambdaL2+LambdaR2)/(24.0*pi)
1478 DSIGMADXDY(2)=BBB*CCCv0d*AAA/(AAA**2+
1479 & (Ms2*GGG)**2)+
1480 & BBB*DDDv0d/(Q2-X*Srad-Ms2)
1481 DSIGMADXDY(4)=EEEv0d*
1482 & ((LambdaL2**2+LambdaR2**2)*(1.-Y)**2+
1483 & 2.0*LambdaL2*LambdaR2*Y**2)/
1484 & (AAA**2+Ms2*GAM**2/(24.0*pi)**2)
1485 DSIGMADXDY(3)=FFFv0d*(LambdaL2**2+LambdaR2**2+
1486 & 2.0*Y**2*LambdaL2*LambdaR2)/
1487 & (Q2-X*Srad-Ms2)**2
1488
1489 if(DSIGMADXDY(3)+DSIGMADXDY(4).eq.0)then
1490 goto 999
1491 endif
1492 sfrac1=DSIGMADXDY(4)/(DSIGMADXDY(3)+DSIGMADXDY(4))
1493
1494 rand(1)=pyr(0)
1495 if(rand(1).lt.sfrac1)then
1496
1497 genproc=1
1498 genprtyp=2
1499 if(q_i.eq.1)then
1500 if(echar.gt.0)pvalence=pvalences_d
1501 endif
1502 if(q_j.eq.1)q_o=1
1503 if(q_j.eq.2)q_o=3
1504 if(q_j.eq.3)q_o=5
1505 if(q_i.eq.1)q_s=1
1506 if(q_i.eq.2)q_s=3
1507 if(q_i.eq.3)q_s=5
1508 else
1509
1510 genproc=2
1511 genprtyp=8
1512 if(q_j.eq.1)then
1513 if(echar.lt.0)pvalence=pvalenceu_d
1514 endif
1515 if(q_i.eq.1)q_o=-1
1516 if(q_i.eq.2)q_o=-3
1517 if(q_i.eq.3)q_o=-5
1518 if(q_j.eq.1)q_s=-1
1519 if(q_j.eq.2)q_s=-3
1520 if(q_j.eq.3)q_s=-5
1521 endif
1522 ELSE IF (LQTYPE.EQ.10) THEN
1523 AAA=(X*Srad-Ms2)
1524 GAM=MLQ*(sqrt(2.0)*LambdaL2_1+LambdaR2_1)**2
1525 if(l_o.gt.1)GAM=GAM+
1526 & MLQ*(sqrt(2.0)*LambdaL2_2+LambdaR2_2)**2
1527 BBB=-LambdaR2*B_RR_U-2.0*LambdaL2*B_LL_U
1528 GGG=(2.0*LambdaL2+LambdaR2)/(24.0*pi)
1529 DSIGMADXDY(2)=BBB*CCCv0u*AAA/(AAA**2+
1530 & (Ms2*GGG)**2)+
1531 & BBB*DDDv0u/(Q2-X*Srad-Ms2)
1532 if(q_j.ne.3)then
1533 DSIGMADXDY(4)=EEEv0u*
1534 & ((4.0*LambdaL2**2+LambdaR2**2)*(1.-Y)**2+
1535 & 2.0*2.0*LambdaL2*LambdaR2*Y**2)/
1536 & (AAA**2+Ms2*GAM**2/(24.0*pi)**2)
1537 else
1538
1539 DSIGMADXDY(4)=0.d0
1540 endif
1541 if(q_i.ne.3)then
1542 DSIGMADXDY(3)=
1543 & FFFv0u*(4.0*LambdaL2**2+LambdaR2**2+
1544 & 2.0*Y**2*2.0*LambdaL2*LambdaR2)/
1545 & (Q2-X*Srad-Ms2)**2
1546 else
1547
1548 DSIGMADXDY(3)=0.d0
1549 endif
1550
1551 if(DSIGMADXDY(3)+DSIGMADXDY(4).eq.0)then
1552 goto 999
1553 endif
1554 sfrac1=DSIGMADXDY(4)/(DSIGMADXDY(3)+DSIGMADXDY(4))
1555
1556 rand(1)=pyr(0)
1557 if(rand(1).lt.sfrac1)then
1558
1559 genproc=1
1560 genprtyp=1
1561 if(q_i.eq.1)then
1562 if(echar.gt.0)pvalence=pvalences_u
1563 endif
1564 if(q_j.eq.1)q_o=2
1565 if(q_j.eq.2)q_o=4
1566 if(q_j.eq.3)q_o=6
1567 if(q_i.eq.1)q_s=2
1568 if(q_i.eq.2)q_s=4
1569 if(q_i.eq.3)q_s=6
1570 else
1571
1572 genproc=2
1573 genprtyp=7
1574 if(q_j.eq.1)then
1575 if(echar.lt.0)pvalence=pvalenceu_u
1576 endif
1577 if(q_i.eq.1)q_o=-2
1578 if(q_i.eq.2)q_o=-4
1579 if(q_i.eq.3)q_o=-6
1580 if(q_j.eq.1)q_s=-2
1581 if(q_j.eq.2)q_s=-4
1582 if(q_j.eq.3)q_s=-6
1583 endif
1584 ELSE IF (LQTYPE.EQ.11) THEN
1585 AAA=(X*Srad-Ms2)
1586 GAM=MLQ*(sqrt(2.0)*LambdaL2_1+LambdaR2_1)**2
1587 if(l_o.gt.1)GAM=GAM+
1588 & MLQ*(sqrt(2.0)*LambdaL2_2+LambdaR2_2)**2
1589 BBB=-LambdaR2*B_RR_D-LambdaL2*B_LL_D
1590 GGG=(2.0*LambdaL2+LambdaR2)/(24.0*pi)
1591 DSIGMADXDY(2)=BBB*CCCv0d*AAA/(AAA**2+
1592 & (Ms2*GGG)**2)+
1593 & BBB*DDDv0d/(Q2-X*Srad-Ms2)
1594 DSIGMADXDY(4)=EEEv0d*
1595 & ((LambdaL2**2+LambdaR2**2)*(1.-Y)**2+
1596 & 2.0*LambdaL2*LambdaR2*Y**2)/
1597 & (AAA**2+Ms2*GAM**2/(24.0*pi)**2)
1598 DSIGMADXDY(3)=FFFv0d*(LambdaL2**2+LambdaR2**2+
1599 & 2.0*Y**2*LambdaL2*LambdaR2)/
1600 & (Q2-X*Srad-Ms2)**2
1601
1602 sfrac1=DSIGMADXDY(4)
1603 ufrac1=DSIGMADXDY(3)
1604 BBB=-LambdaR2*B_RR_U-2.0*LambdaL2*B_LL_U
1605 GGG=(2.0*LambdaL2+LambdaR2)/(24.0*pi)
1606 DSIGMADXDY(2)=DSIGMADXDY(2)+BBB*CCCv0u*AAA/(AAA**2+
1607 & (Ms2*GGG)**2)+
1608 & BBB*DDDv0u/(Q2-X*Srad-Ms2)
1609 if(q_j.ne.3)then
1610
1611 DSIGMADXDY(4)=DSIGMADXDY(4)+EEEv0u*
1612 & ((4.0*LambdaL2**2+LambdaR2**2)*(1.-Y)**2+
1613 & 2.0*2.0*LambdaL2*LambdaR2*Y**2)/
1614 & (AAA**2+Ms2*GAM**2/(24.0*pi)**2)
1615 endif
1616 if(q_i.ne.3)then
1617
1618 DSIGMADXDY(3)=DSIGMADXDY(3)+
1619 & FFFv0u*(4.0*LambdaL2**2+LambdaR2**2+
1620 & 2.0*Y**2*2.0*LambdaL2*LambdaR2)/
1621 & (Q2-X*Srad-Ms2)**2
1622 endif
1623
1624 if(DSIGMADXDY(3)+DSIGMADXDY(4).eq.0)then
1625 goto 999
1626 endif
1627 sfrac2=DSIGMADXDY(4)-sfrac1
1628 ufrac2=DSIGMADXDY(3)-ufrac1
1629 sfrac1=sfrac1/(DSIGMADXDY(3)+DSIGMADXDY(4))
1630 sfrac2=sfrac2/(DSIGMADXDY(3)+DSIGMADXDY(4))
1631 ufrac1=ufrac1/(DSIGMADXDY(3)+DSIGMADXDY(4))
1632 ufrac2=ufrac2/(DSIGMADXDY(3)+DSIGMADXDY(4))
1633
1634 rand(1)=pyr(0)
1635 if(rand(1).lt.sfrac1)then
1636
1637 genproc=1
1638 genprtyp=2
1639 if(q_i.eq.1)then
1640 if(echar.gt.0)pvalence=pvalences_d
1641 endif
1642 if(q_j.eq.1)q_o=1
1643 if(q_j.eq.2)q_o=3
1644 if(q_j.eq.3)q_o=5
1645 if(q_i.eq.1)q_s=1
1646 if(q_i.eq.2)q_s=3
1647 if(q_i.eq.3)q_s=5
1648 elseif(rand(1).lt.(sfrac1+sfrac2))then
1649
1650 genproc=1
1651 genprtyp=1
1652 if(q_i.eq.1)then
1653 if(echar.gt.0)pvalence=pvalences_u
1654 endif
1655 if(q_j.eq.1)q_o=2
1656 if(q_j.eq.2)q_o=4
1657 if(q_j.eq.3)q_o=6
1658 if(q_i.eq.1)q_s=2
1659 if(q_i.eq.2)q_s=4
1660 if(q_i.eq.3)q_s=6
1661 elseif(rand(1).lt.(sfrac1+sfrac2+ufrac1))then
1662
1663 genproc=2
1664 genprtyp=8
1665 if(q_j.eq.1)then
1666 if(echar.lt.0)pvalence=pvalenceu_d
1667 endif
1668 if(q_i.eq.1)q_o=-1
1669 if(q_i.eq.2)q_o=-3
1670 if(q_i.eq.3)q_o=-5
1671 if(q_j.eq.1)q_s=-1
1672 if(q_j.eq.2)q_s=-3
1673 if(q_j.eq.3)q_s=-5
1674 else
1675
1676 genproc=2
1677 genprtyp=7
1678 if(q_j.eq.1)then
1679 if(echar.lt.0)pvalence=pvalenceu_u
1680 endif
1681 if(q_i.eq.1)q_o=-2
1682 if(q_i.eq.2)q_o=-4
1683 if(q_i.eq.3)q_o=-6
1684 if(q_j.eq.1)q_s=-2
1685 if(q_j.eq.2)q_s=-4
1686 if(q_j.eq.3)q_s=-6
1687 endif
1688 ELSE IF (LQTYPE.EQ.12) THEN
1689 GAM=MLQ*(LambdaL2_1+LambdaR2_1)**2
1690 if(l_o.gt.1)GAM=GAM+
1691 & MLQ*(LambdaL2_2+LambdaR2_2)**2
1692 AAA=(X*Srad-Ms2)**2+
1693 & Ms2*GAM**2/((16.0*pi)**2)
1694 BBB=-LambdaR2*B_RL_U-LambdaL2*B_LR_U
1695 DSIGMADXDY(2)=BBB*CCCf0u/(Q2-X*Srad-Ms2)+
1696 & BBB*DDDf0u*(X*Srad-Ms2)/AAA
1697 GGG=(LambdaR2+LambdaL2)**2
1698 if(q_i.ne.3)then
1699 DSIGMADXDY(3)=EEEf0u*GGG/(Q2-X*Srad-Ms2)**2
1700 else
1701
1702 DSIGMADXDY(3)=0.d0
1703 endif
1704 if(q_j.ne.3)then
1705 DSIGMADXDY(4)=FFFf0u*GGG/AAA
1706 else
1707
1708 DSIGMADXDY(4)=0.d0
1709 endif
1710
1711 if(DSIGMADXDY(3)+DSIGMADXDY(4).eq.0)then
1712 goto 999
1713 endif
1714 if(DSIGMADXDY(3).gt.0)then
1715 sfrac1=DSIGMADXDY(4)/(DSIGMADXDY(3)+DSIGMADXDY(4))
1716 else
1717 sfrac1=0.d0
1718 endif
1719
1720 rand(1)=pyr(0)
1721 if(rand(1).lt.sfrac1)then
1722
1723 genproc=1
1724 genprtyp=1
1725 if(q_i.eq.1)then
1726 if(echar.gt.0)pvalence=pvalences_u
1727 endif
1728 if(q_j.eq.1)q_o=2
1729 if(q_j.eq.2)q_o=4
1730 if(q_j.eq.3)q_o=6
1731 if(q_i.eq.1)q_s=2
1732 if(q_i.eq.2)q_s=4
1733 if(q_i.eq.3)q_s=6
1734 else
1735
1736 genproc=2
1737 genprtyp=7
1738 if(q_j.eq.1)then
1739 if(echar.lt.0)pvalence=pvalenceu_u
1740 endif
1741 if(q_i.eq.1)q_o=-2
1742 if(q_i.eq.2)q_o=-4
1743 if(q_i.eq.3)q_o=-6
1744 if(q_j.eq.1)q_s=-2
1745 if(q_j.eq.2)q_s=-4
1746 if(q_j.eq.3)q_s=-6
1747 endif
1748 ELSE IF (LQTYPE.EQ.13) THEN
1749 GAM=MLQ*(LambdaL2_1+LambdaR2_1)**2
1750 if(l_o.gt.1)GAM=GAM+
1751 & MLQ*(LambdaL2_2+LambdaR2_2)**2
1752 AAA=(X*Srad-Ms2)**2+
1753 & Ms2*GAM**2/((16.0*pi)**2)
1754 BBB=-LambdaR2*B_RL_U-LambdaL2*B_LR_U
1755 DSIGMADXDY(2)=BBB*CCCf0u/(Q2-X*Srad-Ms2)+
1756 & BBB*DDDf0u*(X*Srad-Ms2)/AAA
1757 GGG=(LambdaR2+LambdaL2)**2
1758 if(q_i.ne.3)then
1759 DSIGMADXDY(3)=EEEf0u*GGG/(Q2-X*Srad-Ms2)**2
1760 else
1761
1762 DSIGMADXDY(3)=0.d0
1763 endif
1764 if(q_j.ne.3)then
1765 DSIGMADXDY(4)=FFFf0u*GGG/AAA
1766 else
1767
1768 DSIGMADXDY(3)=0.d0
1769 endif
1770
1771 sfrac1=DSIGMADXDY(4)
1772 ufrac1=DSIGMADXDY(3)
1773 BBB=-LambdaR2*B_RL_D-LambdaL2*B_LR_D
1774 DSIGMADXDY(2)=DSIGMADXDY(2)+BBB*CCCf0d/(Q2-X*Srad-Ms2)+
1775 & BBB*DDDf0d*(X*Srad-Ms2)/AAA
1776 GGG=(LambdaR2+LambdaL2)**2
1777 DSIGMADXDY(3)=DSIGMADXDY(3)+EEEf0d*GGG/(Q2-X*Srad-Ms2)**2
1778 DSIGMADXDY(4)=DSIGMADXDY(4)+FFFf0d*GGG/AAA
1779 if(DSIGMADXDY(3)+DSIGMADXDY(4).eq.0)then
1780 goto 999
1781 endif
1782 sfrac2=DSIGMADXDY(4)-sfrac1
1783 ufrac2=DSIGMADXDY(3)-ufrac1
1784 sfrac1=sfrac1/(DSIGMADXDY(3)+DSIGMADXDY(4))
1785 sfrac2=sfrac2/(DSIGMADXDY(3)+DSIGMADXDY(4))
1786 ufrac1=ufrac1/(DSIGMADXDY(3)+DSIGMADXDY(4))
1787 ufrac2=ufrac2/(DSIGMADXDY(3)+DSIGMADXDY(4))
1788
1789 rand(1)=pyr(0)
1790 if(rand(1).lt.sfrac1)then
1791
1792 genproc=1
1793 genprtyp=1
1794 if(q_i.eq.1)then
1795 if(echar.gt.0)pvalence=pvalences_u
1796 endif
1797 if(q_j.eq.1)q_o=2
1798 if(q_j.eq.2)q_o=4
1799 if(q_j.eq.3)q_o=6
1800 if(q_i.eq.1)q_s=2
1801 if(q_i.eq.2)q_s=4
1802 if(q_i.eq.3)q_s=6
1803 elseif(rand(1).lt.(sfrac1+sfrac2))then
1804
1805 genproc=1
1806 genprtyp=2
1807 if(q_i.eq.1)then
1808 if(echar.gt.0)pvalence=pvalences_d
1809 endif
1810 if(q_j.eq.1)q_o=1
1811 if(q_j.eq.2)q_o=3
1812 if(q_j.eq.3)q_o=5
1813 if(q_i.eq.1)q_s=1
1814 if(q_i.eq.2)q_s=3
1815 if(q_i.eq.3)q_s=5
1816 elseif(rand(1).lt.(sfrac1+sfrac2+ufrac1))then
1817
1818 genproc=2
1819 genprtyp=7
1820 if(q_j.eq.1)then
1821 if(echar.lt.0)pvalence=pvalenceu_u
1822 endif
1823 if(q_i.eq.1)q_o=-2
1824 if(q_i.eq.2)q_o=-4
1825 if(q_i.eq.3)q_o=-5
1826 if(q_j.eq.1)q_s=-2
1827 if(q_j.eq.2)q_s=-4
1828 if(q_j.eq.3)q_s=-5
1829 else
1830
1831 genproc=2
1832 genprtyp=8
1833 if(q_j.eq.1)then
1834 if(echar.lt.0)pvalence=pvalenceu_d
1835 endif
1836 if(q_i.eq.1)q_o=-1
1837 if(q_i.eq.2)q_o=-3
1838 if(q_i.eq.3)q_o=-5
1839 if(q_j.eq.1)q_s=-1
1840 if(q_j.eq.2)q_s=-3
1841 if(q_j.eq.3)q_s=-5
1842 endif
1843 ELSE IF (LQTYPE.EQ.14) THEN
1844 GAM=MLQ*(LambdaL2_1+LambdaR2_1)**2
1845 if(l_o.gt.1)GAM=GAM+
1846 & MLQ*(LambdaL2_2+LambdaR2_2)**2
1847 AAA=(X*Srad-Ms2)**2+
1848 & Ms2*GAM**2/((16.0*pi)**2)
1849 BBB=-LambdaR2*B_RL_D-LambdaL2*B_LR_D
1850 DSIGMADXDY(2)=BBB*CCCf0d/(Q2-X*Srad-Ms2)+
1851 & BBB*DDDf0d*(X*Srad-Ms2)/AAA
1852 GGG=(LambdaR2+LambdaL2)**2
1853 DSIGMADXDY(3)=EEEf0d*GGG/(Q2-X*Srad-Ms2)**2
1854 DSIGMADXDY(4)=FFFf0d*GGG/AAA
1855
1856 if(DSIGMADXDY(3)+DSIGMADXDY(4).eq.0)then
1857 goto 999
1858 endif
1859 sfrac1=DSIGMADXDY(4)/(DSIGMADXDY(3)+DSIGMADXDY(4))
1860
1861 rand(1)=pyr(0)
1862 if(rand(1).lt.sfrac1)then
1863
1864 genproc=1
1865 genprtyp=2
1866 if(q_i.eq.1)then
1867 if(echar.gt.0)pvalence=pvalences_d
1868 endif
1869 if(q_j.eq.1)q_o=1
1870 if(q_j.eq.2)q_o=3
1871 if(q_j.eq.3)q_o=5
1872 if(q_i.eq.1)q_s=1
1873 if(q_i.eq.2)q_s=3
1874 if(q_i.eq.3)q_s=5
1875 else
1876
1877 genproc=2
1878 genprtyp=8
1879 if(q_j.eq.1)then
1880 if(echar.lt.0)pvalence=pvalenceu_d
1881 endif
1882 if(q_i.eq.1)q_o=-1
1883 if(q_i.eq.2)q_o=-3
1884 if(q_i.eq.3)q_o=-5
1885 if(q_j.eq.1)q_s=-1
1886 if(q_j.eq.2)q_s=-3
1887 if(q_j.eq.3)q_s=-5
1888 endif
1889 ENDIF
1890 DSIGMADXDY(1)=DSIGMADXDY(1)*pi*alpha**2/(Srad*(X*Y)**2)
1891 DSIGMADXDY(2)=DSIGMADXDY(2)*pi*alpha**2/(Srad*(X*Y)**2)
1892 DSIGMADXDY(3)=DSIGMADXDY(3)*pi*alpha**2/(Srad*(X*Y)**2)
1893 DSIGMADXDY(4)=DSIGMADXDY(4)*pi*alpha**2/(Srad*(X*Y)**2)
1894 DSIGMADXDY(1)=DSIGMADXDY(1)*0.38938*(10.0)**9
1895 DSIGMADXDY(2)=DSIGMADXDY(2)*0.38938*(10.0)**9
1896 DSIGMADXDY(3)=DSIGMADXDY(3)*0.38938*(10.0)**9
1897 DSIGMADXDY(4)=DSIGMADXDY(4)*0.38938*(10.0)**9
1898
1899 999 continue
1900
1901 if(echar.lt.0)q_o=-q_o
1902 if(echar.lt.0)q_s=-q_s
1903 if(l_o.eq.1)then
1904
1905 DXSEC(1)=DSIGMADXDY(2)
1906 else
1907 DXSEC(1)=0.d0
1908 endif
1909 DXSEC(2)=DSIGMADXDY(3)
1910 DXSEC(3)=DSIGMADXDY(4)
1911
1912 RETURN
1913 END
1914
1915 subroutine LQGBAN
1916
1917
1918
1919 implicit none
1920
1921 write(6,*)
1922 write(6,*) ' *************************************************'
1923 write(6,*) '- -'
1924 write(6,*) '- LQGENEP -'
1925 write(6,*) '- -'
1926 write(6,*) '- -'
1927 write(6,*) '- -'
1928 write(6,*) '- LeptoQuark GENerator for Electron-Proton -'
1929 write(6,*) '- scattering -'
1930 write(6,*) '- -'
1931 write(6,*) '- -'
1932 write(6,*) '- -'
1933 write(6,*) '- Author: L.Bellagamba -'
1934 write(6,*) '- e-mail lorenzo.bellagamba@bo.infn.it -'
1935 write(6,*) '- -'
1936 write(6,*) '- Version: 1.0 -'
1937 write(6,*) '- Date: 01.03.2001 -'
1938 write(6,*) ' *************************************************'
1939 write(6,*)
1940 write(6,*)
1941 write(6,*)
1942 return
1943 end
1944
1945 subroutine LQGPRI1
1946
1947
1948
1949
1950 implicit none
1951
1952
1953 double precision BEAMPAR,LQGPAR3
1954 integer LQGPAR1,LQGPAR2
1955 CHARACTER*256 out_file
1956 COMMON/LQGDATA/BEAMPAR(3),LQGPAR1(10),LQGPAR2(10),LQGPAR3(20),
1957 > out_file
1958
1959
1960 write(6,*) '>>>>>>>>>>>>>> LQGENEP RUN REQUEST <<<<<<<<<<<<<<'
1961 write(6,*)
1962 write(6,101) LQGPAR1(1)
1963 write(6,1011) BEAMPAR(1)
1964 write(6,1012) BEAMPAR(2)
1965 write(6,1013) BEAMPAR(3)
1966 write(6,102) LQGPAR2(1)
1967 write(6,103) sngl(LQGPAR3(1))
1968 write(6,104) sngl(LQGPAR3(2))
1969 write(6,105) sngl(LQGPAR3(3))
1970 write(6,106) LQGPAR2(2)
1971 write(6,107) LQGPAR2(3)
1972 write(6,108) LQGPAR2(4)
1973 write(6,*)
1974 write(6,109) sngl(LQGPAR3(4)),sngl(LQGPAR3(5))
1975 write(6,110) sngl(LQGPAR3(6)),sngl(LQGPAR3(7))
1976 write(6,111) sngl(LQGPAR3(8))
1977 write(6,*)
1978 write(6,*) '>>>>>>>>>>>>>>>>>>>>>>>>> <<<<<<<<<<<<<<<<<<<<<<<<<'
1979 101 format(' Number of events.................. ',I8)
1980 1011 format(' Electron beam charge.............. ',F3.0)
1981 1012 format(' Electron beam energy.............. ',F7.2)
1982 1013 format(' Proton beam energy................ ',F7.2)
1983 102 format(' LQ type........................... ',I2)
1984 103 format(' LQ mass (GeV)..................... ',F7.2)
1985 104 format(' LQ production coupling (s-ch.).... ',F7.4)
1986 105 format(' LQ decay coupling (s-ch.)......... ',F7.4)
1987 106 format(' struck quark generation (s-ch.)... ',I2)
1988 107 format(' output quark generation (s-ch.)... ',I2)
1989 108 format(' output lepton generation.......... ',I2)
1990 109 format(' X generation range................ ',F6.3,' - ',F6.3)
1991 110 format(' Y generation range................ ',F6.3,' - ',F6.3)
1992 111 format(' minimum allowed Q2 (GeV^2)........ ',F7.2)
1993 return
1994 end
1995
1996 subroutine LQGPRI2
1997
1998
1999
2000
2001 Implicit none
2002
2003 double precision DXSEC(3),pvalence
2004 integer q_o,q_s,genproc,genprtyp,sch,uch,gproc(8)
2005 common /LQGout/ DXSEC,pvalence,q_o,q_s,genproc,genprtyp,
2006 >sch,uch,gproc
2007
2008
2009 double precision BEAMPAR,LQGPAR3
2010 integer LQGPAR1,LQGPAR2
2011 CHARACTER*256 out_file
2012 COMMON/LQGDATA/BEAMPAR(3),LQGPAR1(10),LQGPAR2(10),LQGPAR3(20),
2013 > out_file
2014
2015
2016
2017 double precision LQGKPAR,LQGDDX
2018 integer LQGPID
2019 COMMON/LQGEVT/LQGKPAR(3),LQGDDX(3),LQGPID(3)
2020
2021
2022 double precision Mlq,G1,G2
2023 Integer LQtype,l_o,q_i,q_j
2024 common /LQGproc/ Mlq,G1,G2,LQtype,l_o,q_i,q_j
2025
2026
2027 integer echar
2028 double precision ebeam,pbeam
2029 common /LQGbeam/ebeam,pbeam,echar
2030
2031 Character*40 CHANS1,CHANS2,CHANU1,CHANU2
2032
2033 CALL LQGCHAN(CHANS1,CHANS2,CHANU1,CHANU2)
2034
2035 write(6,*) '>>>>>>>>>>> LQGENEP FINAL STATISTICS <<<<<<<<<<<'
2036 write(6,*)
2037 write(6,101) LQGPAR1(4)
2038 write(6,104) sch
2039 write(6,106) uch
2040 write(6,*)
2041 write(6,*)
2042 write(6,*) ' Details of the generation '
2043 write(6,*) ' --------------------------- '
2044 write(6,*)
2045 write(6,*) '=================================================='
2046 write(6,*) 'I I I'
2047 write(6,*) 'I process I events I'
2048 write(6,*) 'I I I'
2049 write(6,*) '=================================================='
2050 write(6,*) 'I I I'
2051 if(lqtype.eq.1.or.lqtype.eq.2)then
2052 write(6,111)
2053 write(6,110) CHANS1,gproc(3)
2054 write(6,*) 'I ----------------------- I I'
2055 write(6,112)
2056 write(6,110) CHANU1,gproc(5)
2057 elseif(lqtype.eq.3)then
2058 write(6,111)
2059 write(6,110) CHANS1,gproc(4)
2060 write(6,*) 'I ----------------------- I I'
2061 write(6,112)
2062 write(6,110) CHANU1,gproc(6)
2063 elseif(lqtype.eq.4)then
2064 write(6,111)
2065 write(6,110) CHANS1,gproc(4)
2066 write(6,110) CHANS2,gproc(3)
2067 write(6,*) 'I ----------------------- I I'
2068 write(6,112)
2069 write(6,110) CHANU1,gproc(6)
2070 write(6,110) CHANU2,gproc(5)
2071 elseif(lqtype.eq.5)then
2072 write(6,111)
2073 write(6,110) CHANS1,gproc(4)
2074 write(6,*) 'I ----------------------- I I'
2075 write(6,112)
2076 write(6,110) CHANU1,gproc(6)
2077 elseif(lqtype.eq.6)then
2078 write(6,111)
2079 write(6,110) CHANS1,gproc(4)
2080 write(6,110) CHANS2,gproc(3)
2081 write(6,*) 'I ----------------------- I I'
2082 write(6,112)
2083 write(6,110) CHANU1,gproc(6)
2084 write(6,110) CHANU2,gproc(5)
2085 elseif(lqtype.eq.7)then
2086 write(6,111)
2087 write(6,110) CHANS1,gproc(3)
2088 write(6,*) 'I ----------------------- I I'
2089 write(6,112)
2090 write(6,110) CHANU1,gproc(5)
2091 elseif(lqtype.eq.8.or.lqtype.eq.9)then
2092 write(6,111)
2093 write(6,110) CHANS1,gproc(2)
2094 write(6,*) 'I ----------------------- I I'
2095 write(6,112)
2096 write(6,110) CHANU1,gproc(8)
2097 elseif(lqtype.eq.10)then
2098 write(6,111)
2099 write(6,110) CHANS1,gproc(1)
2100 write(6,*) 'I ----------------------- I I'
2101 write(6,112)
2102 write(6,110) CHANU1,gproc(7)
2103 elseif(lqtype.eq.11)then
2104 write(6,111)
2105 write(6,110) CHANS1,gproc(2)
2106 write(6,110) CHANS2,gproc(1)
2107 write(6,*) 'I ----------------------- I I'
2108 write(6,112)
2109 write(6,110) CHANU1,gproc(8)
2110 write(6,110) CHANU2,gproc(7)
2111 elseif(lqtype.eq.12)then
2112 write(6,111)
2113 write(6,110) CHANS1,gproc(1)
2114 write(6,*) 'I ----------------------- I I'
2115 write(6,112)
2116 write(6,110) CHANU1,gproc(7)
2117 elseif(lqtype.eq.13)then
2118 write(6,111)
2119 write(6,110) CHANS1,gproc(2)
2120 write(6,110) CHANS2,gproc(1)
2121 write(6,*) 'I ----------------------- I I'
2122 write(6,112)
2123 write(6,110) CHANU1,gproc(8)
2124 write(6,110) CHANU2,gproc(7)
2125 elseif(lqtype.eq.14)then
2126 write(6,111)
2127 write(6,110) CHANS1,gproc(2)
2128 write(6,*) 'I ----------------------- I I'
2129 write(6,112)
2130 write(6,110) CHANU1,gproc(8)
2131 endif
2132 write(6,*) 'I I I'
2133 write(6,*) '=================================================='
2134 write(6,*)
2135 write(6,*) ' --------------------------- '
2136 write(6,*)
2137 write(6,*)
2138 write(6,*) '>>>>>>>>>>>>>>>>>>>>>>>>> <<<<<<<<<<<<<<<<<<<<<<<<'
2139 101 format(' Number of generated events.............',I8)
2140
2141
2142 104 format(' Number of s-channel generated events...',I8)
2143
2144 106 format(' Number of u-channel generated events...',I8)
2145 110 format(1X,'I ',A38,'I',1X,I6,' I')
2146 111 format(1X,'I ','s-channel:',28X,'I',7X,' I')
2147 112 format(1X,'I ','u-channel:',28X,'I',7X,' I')
2148 return
2149 end
2150
2151 Subroutine LQGCHAN(CHANS1,CHANS2,CHANU1,CHANU2)
2152
2153
2154
2155
2156
2157 Implicit None
2158
2159 Character*38 CHANS1,CHANS2,CHANU1,CHANU2
2160
2161
2162 integer echar
2163 double precision ebeam,pbeam
2164 common /LQGbeam/ebeam,pbeam,echar
2165
2166
2167 double precision Mlq,G1,G2
2168 Integer LQtype,l_o,q_i,q_j
2169 common /LQGproc/ Mlq,G1,G2,LQtype,l_o,q_i,q_j
2170
2171 Character*7 LQCHA(14)
2172 DATA LQCHA /'S_0L','S_0R','~S_0R','S_1L',
2173 > 'V_1/2L','V_1/2R','~V_1/2L',
2174 > 'V_0L','V_0R','~V_0R','V_1L',
2175 > 'S_1/2L','S_1/2R','~S_1/2L'/
2176 character*5 uch(3), dch(3)
2177 character*5 ubch(3), dbch(3)
2178 character*5 lchp(3),lchm(3)
2179 data uch,ubch /' u ',' c ',' t ', ' ubar',' cbar',' tbar'/
2180 data dch,dbch /' d ',' s ',' b ', ' dbar',' sbar',' bbar'/
2181 data lchp,lchm /'e+','mu+','tau+','e-','mu-','tau-'/
2182 character*5 l_in
2183 character*5 q_in
2184 character*5 l_ou
2185 character*5 q_ou
2186
2187 if(echar.gt.0) then
2188 l_in ='e+'
2189 l_ou=lchp(l_o)
2190 else
2191 l_in ='e-'
2192 l_ou=lchm(l_o)
2193 endif
2194 IF(lqtype.eq.1.or.lqtype.eq.2.or.lqtype.eq.7)Then
2195 if(echar.gt.0) then
2196 q_in=ubch(q_i)
2197 q_ou=ubch(q_j)
2198 else
2199 q_in=uch(q_i)
2200 q_ou=uch(q_j)
2201 endif
2202 CHANS1=l_in//q_in//' -> '//LQCHA(lqtype)//' -> '//l_ou//q_ou
2203 CHANS2=' '
2204 if(echar.gt.0) then
2205 q_in=uch(q_j)
2206 q_ou=uch(q_i)
2207 else
2208 q_in=ubch(q_i)
2209 q_ou=ubch(q_j)
2210 endif
2211 CHANU1=l_in//q_in//' -> '//LQCHA(lqtype)//' -> '//l_ou//q_ou
2212 CHANU2=' '
2213 ELSEIF(lqtype.eq.3.or.lqtype.eq.5)Then
2214 if(echar.gt.0) then
2215 q_in=dbch(q_i)
2216 q_ou=dbch(q_j)
2217 else
2218 q_in=dch(q_i)
2219 q_ou=dch(q_j)
2220 endif
2221 CHANS1=l_in//q_in//' -> '//LQCHA(lqtype)//' -> '//l_ou//q_ou
2222 CHANS2=' '
2223 if(echar.gt.0) then
2224 q_in=dch(q_j)
2225 q_ou=dch(q_i)
2226 else
2227 q_in=dbch(q_j)
2228 q_ou=dbch(q_i)
2229 endif
2230 CHANU1=l_in//q_in//' -> '//LQCHA(lqtype)//' -> '//l_ou//q_ou
2231 CHANU2=' '
2232 ELSEIF(lqtype.eq.4.or.lqtype.eq.6)Then
2233 if(echar.gt.0) then
2234 q_in=dbch(q_i)
2235 q_ou=dbch(q_j)
2236 else
2237 q_in=dch(q_i)
2238 q_ou=dch(q_j)
2239 endif
2240 CHANS1=l_in//q_in//' -> '//LQCHA(lqtype)//' -> '//l_ou//q_ou
2241 if(echar.gt.0) then
2242 q_in=ubch(q_i)
2243 q_ou=ubch(q_j)
2244 else
2245 q_in=uch(q_i)
2246 q_ou=uch(q_j)
2247 endif
2248 CHANS2=l_in//q_in//' -> '//LQCHA(lqtype)//' -> '//l_ou//q_ou
2249 if(echar.gt.0) then
2250 q_in=dch(q_j)
2251 q_ou=dch(q_i)
2252 else
2253 q_in=dbch(q_j)
2254 q_ou=dbch(q_i)
2255 endif
2256 CHANU1=l_in//q_in//' -> '//LQCHA(lqtype)//' -> '//l_ou//q_ou
2257 if(echar.gt.0) then
2258 q_in=uch(q_j)
2259 q_ou=uch(q_i)
2260 else
2261 q_in=ubch(q_j)
2262 q_ou=ubch(q_i)
2263 endif
2264 CHANU2=l_in//q_in//' -> '//LQCHA(lqtype)//' -> '//l_ou//q_ou
2265 ELSEIF(lqtype.eq.8.or.lqtype.eq.9.or.lqtype.eq.14)Then
2266 if(echar.gt.0) then
2267 q_in=dch(q_i)
2268 q_ou=dch(q_j)
2269 else
2270 q_in=dbch(q_i)
2271 q_ou=dbch(q_j)
2272 endif
2273 CHANS1=l_in//q_in//' -> '//LQCHA(lqtype)//' -> '//l_ou//q_ou
2274 CHANS2=' '
2275 if(echar.gt.0) then
2276 q_in=dbch(q_j)
2277 q_ou=dbch(q_i)
2278 else
2279 q_in=dch(q_i)
2280 q_ou=dch(q_j)
2281 endif
2282 CHANU1=l_in//q_in//' -> '//LQCHA(lqtype)//' -> '//l_ou//q_ou
2283 CHANU2=' '
2284 ELSEIF(lqtype.eq.10.or.lqtype.eq.12)Then
2285 if(echar.gt.0) then
2286 q_in=uch(q_i)
2287 q_ou=uch(q_j)
2288 else
2289 q_in=ubch(q_i)
2290 q_ou=ubch(q_j)
2291 endif
2292 CHANS1=l_in//q_in//' -> '//LQCHA(lqtype)//' -> '//l_ou//q_ou
2293 CHANS2=' '
2294 if(echar.gt.0) then
2295 q_in=ubch(q_j)
2296 q_ou=ubch(q_i)
2297 else
2298 q_in=uch(q_i)
2299 q_ou=uch(q_j)
2300 endif
2301 CHANU1=l_in//q_in//' -> '//LQCHA(lqtype)//' -> '//l_ou//q_ou
2302 CHANU2=' '
2303 ELSEIF(lqtype.eq.11.or.lqtype.eq.13)Then
2304 if(echar.gt.0) then
2305 q_in=dch(q_i)
2306 q_ou=dch(q_j)
2307 else
2308 q_in=dbch(q_i)
2309 q_ou=dbch(q_j)
2310 endif
2311 CHANS1=l_in//q_in//' -> '//LQCHA(lqtype)//' -> '//l_ou//q_ou
2312 if(echar.gt.0) then
2313 q_in=uch(q_i)
2314 q_ou=uch(q_j)
2315 else
2316 q_in=ubch(q_i)
2317 q_ou=ubch(q_j)
2318 endif
2319 CHANS2=l_in//q_in//' -> '//LQCHA(lqtype)//' -> '//l_ou//q_ou
2320 if(echar.gt.0) then
2321 q_in=dbch(q_j)
2322 q_ou=dbch(q_i)
2323 else
2324 q_in=dch(q_j)
2325 q_ou=dch(q_i)
2326 endif
2327 CHANU1=l_in//q_in//' -> '//LQCHA(lqtype)//' -> '//l_ou//q_ou
2328 if(echar.gt.0) then
2329 q_in=ubch(q_j)
2330 q_ou=ubch(q_i)
2331 else
2332 q_in=uch(q_j)
2333 q_ou=uch(q_i)
2334 endif
2335 CHANU2=l_in//q_in//' -> '//LQCHA(lqtype)//' -> '//l_ou//q_ou
2336 ENDIF
2337 Return
2338 End
2339
2340
2341
2342
2343
2344
2345
2346
2347
2348
2349
2350
2351
2352
2353
2354
2355
2356
2357
2358
2359
2360
2361
2362
2363
2364
2365
2366
2367
2368
2369
2370
2371
2372
2373
2374
2375
2376
2377
2378
2379
2380
2381
2382
2383
2384
2385
2386
2387
2388
2389
2390
2391
2392
2393
2394
2395
2396
2397
2398
2399
2400
2401
2402
2403
2404
2405
2406
2407
2408
2409
2410
2411
2412
2413
2414
2415
2416
2417
2418 Subroutine LQGTESTH
2419
2420
2421
2422
2423
2424
2425
2426 implicit none
2427
2428
2429 double precision PARP,PARI
2430 integer MSTP,MSTI
2431 COMMON/PYPARS/MSTP(200),PARP(200),MSTI(200),PARI(200)
2432
2433
2434 double precision BEAMPAR,LQGPAR3
2435 integer LQGPAR1,LQGPAR2
2436 CHARACTER*256 out_file
2437 COMMON/LQGDATA/BEAMPAR(3),LQGPAR1(10),LQGPAR2(10),LQGPAR3(20),
2438 > out_file
2439
2440
2441
2442
2443 double precision LQGKPAR,LQGDDX
2444 integer LQGPID
2445 COMMON/LQGEVT/LQGKPAR(3),LQGDDX(3),LQGPID(3)
2446
2447 integer NeV,i
2448
2449 beampar(1)=1.
2450 beampar(2)=27.5
2451 beampar(3)=820.
2452 lqgpar1(1)=5000
2453 lqgpar1(2)=10
2454 lqgpar1(3)=1
2455 lqgpar1(4)=0
2456 lqgpar1(5)=0
2457 lqgpar2(1)=3
2458 lqgpar2(2)=2
2459 lqgpar2(3)=1
2460 lqgpar2(4)=3
2461 lqgpar3(1)=600.
2462 lqgpar3(2)=0.3
2463 lqgpar3(3)=0.3
2464 lqgpar3(4)=0.
2465 lqgpar3(5)=1.
2466 lqgpar3(6)=0.
2467 lqgpar3(7)=1.
2468 lqgpar3(8)=500.
2469 lqgpar3(9)=1.
2470 lqgpar3(10)=4.
2471 lqgpar3(11)=32.
2472
2473 lqgpar3(12)=2.d-11
2474
2475
2476
2477
2478
2479
2480
2481
2482
2483
2484
2485
2486
2487
2488
2489
2490
2491
2492
2493
2494
2495
2496 return
2497 end
2498
2499 Subroutine LQGXNWA(XNWA,IERR)
2500
2501
2502
2503
2504
2505
2506
2507
2508 Implicit None
2509
2510 double precision XNWA
2511 integer IERR
2512 double precision xx,cutfac0,cutfac1,ycut,factor
2513 double precision upv,dnv,usea,dsea,str,chm,bot,top,gl
2514 double precision xsu,xsd,xsub,xsdb
2515 double precision Br_rat
2516
2517
2518 integer echar
2519 double precision ebeam,pbeam
2520 common /LQGbeam/ebeam,pbeam,echar
2521
2522
2523 double precision xmax,xmin,ymax,ymin,zmax,zmin,Q2min
2524 common /LQGKinC/ xmax,xmin,ymax,ymin,zmax,zmin,Q2min
2525
2526
2527 double precision Mlq,G1,G2
2528 Integer LQtype,l_o,q_i,q_j
2529 common /LQGproc/ Mlq,G1,G2,LQtype,l_o,q_i,q_j
2530
2531
2532 double precision S,Srad,x,y,z,Q2
2533 common /LQGKinV/ S,Srad,x,y,z,Q2
2534
2535 double precision pi
2536 parameter (pi=3.141592653589793d0)
2537
2538 IERR=0
2539 XNWA=0.d0
2540 xx=MLQ**2/S
2541 if (xx.gt.1d0.or.xx.lt.1d-6) then
2542 write(6,*)'#################################################'
2543 write(6,*)'LQGXNWA Error:'
2544 write(6,*)'LQ mass not compatible with center of mass energy'
2545 write(6,*)'#################################################'
2546 IERR=1
2547 return
2548 endif
2549 if (Q2min.le.1.) then
2550 cutfac0=1.d0
2551 cutfac1=2.d0
2552 else
2553 ycut=Q2min/xx/s
2554 if (ycut.ge.1d0) then
2555 XNWA=0.d0
2556 else
2557 cutfac0=1.d0-ycut
2558 cutfac1=2.*(cutfac0)**3
2559 endif
2560 endif
2561 factor=(pi*G1**2)/(4.d0*s)*(0.3894)
2562 call structm(xx,mlq,upv,dnv,usea,dsea,str,chm,bot,top,gl)
2563 if (echar.ge.0) then
2564 xsu=factor*(upv+usea)/xx
2565 xsd=factor*(dnv+dsea)/xx
2566 xsub=factor*usea/xx
2567 xsdb=factor*dsea/xx
2568 else
2569 xsu=factor*usea/xx
2570 xsd=factor*dsea/xx
2571 xsub=factor*(upv+usea)/xx
2572 xsdb=factor*(dnv+dsea)/xx
2573 endif
2574 IF(l_o.eq.1)then
2575
2576 IF(lqtype.eq.1)XNWA = cutfac0*xsub*0.5
2577 IF(lqtype.eq.2)XNWA = cutfac0*xsub
2578 IF(lqtype.eq.3)XNWA = cutfac0*xsdb
2579 IF(lqtype.eq.4)XNWA = cutfac0*(xsub*0.5+2.*xsdb)
2580 IF(lqtype.eq.12)XNWA = cutfac0*xsu
2581 IF(lqtype.eq.13)XNWA = cutfac0*(xsu+xsd)
2582 IF(lqtype.eq.14)XNWA = cutfac0*xsd
2583 IF(lqtype.eq.5)XNWA = cutfac1*(xsdb)
2584 IF(lqtype.eq.6)XNWA = cutfac1*(xsub+xsdb)
2585 IF(lqtype.eq.7)XNWA = cutfac1*(xsub)*0.5
2586 IF(lqtype.eq.8)XNWA = cutfac1*(xsd)*0.5
2587 IF(lqtype.eq.9)XNWA = cutfac1*(xsd)
2588 IF(lqtype.eq.10)XNWA = cutfac1*(xsu)
2589 IF(lqtype.eq.11)XNWA = cutfac1*(2.*xsu+xsd*0.5)
2590 ELSE
2591
2592 Br_rat=G2**2/(G1**2+G2**2)
2593 IF(lqtype.eq.1)XNWA = Br_rat*cutfac0*xsub*0.5
2594 IF(lqtype.eq.2)XNWA = Br_rat*cutfac0*xsub
2595 IF(lqtype.eq.3)XNWA = Br_rat*cutfac0*xsdb
2596 IF(lqtype.eq.4)XNWA = Br_rat*cutfac0*(xsub*0.5+2.*xsdb)
2597 IF(lqtype.eq.12)XNWA = Br_rat*cutfac0*xsu
2598 IF(lqtype.eq.13)XNWA = Br_rat*cutfac0*(xsu+xsd)
2599 IF(lqtype.eq.14)XNWA = Br_rat*cutfac0*xsd
2600 IF(lqtype.eq.5)XNWA = Br_rat*cutfac1*(xsdb)
2601 IF(lqtype.eq.6)XNWA = Br_rat*cutfac1*(xsub+xsdb)
2602 IF(lqtype.eq.7)XNWA = Br_rat*cutfac1*(xsub)*0.5
2603 IF(lqtype.eq.8)XNWA = Br_rat*cutfac1*(xsd)*0.5
2604 IF(lqtype.eq.9)XNWA = Br_rat*cutfac1*(xsd)
2605 IF(lqtype.eq.10)XNWA = Br_rat*cutfac1*(xsu)
2606 IF(lqtype.eq.11)XNWA = Br_rat*cutfac1*(2.*xsu+xsd*0.5)
2607 ENDIF
2608 return
2609 end
2610
2611 subroutine LQGXINT(XINT,XINTE,IERR)
2612
2613
2614
2615
2616
2617
2618
2619
2620
2621 Implicit None
2622
2623 Double precision XINT,XINTE
2624 Integer IERR
2625
2626 double precision relerr,result
2627 external DSDXDY
2628 integer minpts,maxpts,iwk,nfnevl,ifail
2629 parameter (minpts=10000)
2630 parameter (maxpts=5000000)
2631 parameter (iwk=1000000)
2632 double precision wk(iwk),a(2),b(2),eps
2633 parameter (eps=1d-3)
2634
2635
2636 double precision xmax,xmin,ymax,ymin,zmax,zmin,Q2min
2637 common /LQGKinC/ xmax,xmin,ymax,ymin,zmax,zmin,Q2min
2638
2639 IERR=0
2640 XINT=0.d0
2641 XINTE=0.d0
2642
2643 a(1)=XMIN
2644 b(1)=XMAX
2645 a(2)=YMIN
2646 b(2)=YMAX
2647
2648 call dadmul(DSDXDY,2,a,b,minpts,maxpts,eps,
2649 + wk,iwk,result,relerr,nfnevl,ifail)
2650
2651 if(ifail.ne.0)then
2652 Write(6,*)'LQGXINT: Error in Cross Section integration'
2653 IERR=ifail
2654 return
2655 endif
2656
2657 XINT=result
2658 XINTE=abs(result)*relerr
2659
2660 return
2661 end
2662
2663 Double precision function dsdxdy(n,xx)
2664
2665 IMPLICIT None
2666
2667 integer n
2668 double precision xx(*)
2669
2670
2671 double precision xmax,xmin,ymax,ymin,zmax,zmin,Q2min
2672 common /LQGKinC/ xmax,xmin,ymax,ymin,zmax,zmin,Q2min
2673
2674
2675 double precision Mlq,G1,G2
2676 Integer LQtype,l_o,q_i,q_j
2677 common /LQGproc/ Mlq,G1,G2,LQtype,l_o,q_i,q_j
2678
2679
2680 double precision S,Srad,x,y,z,Q2
2681 common /LQGKinV/ S,Srad,x,y,z,Q2
2682
2683
2684 double precision DXSEC(3),pvalence
2685 integer q_o,q_s,genproc,genprtyp,sch,uch,gproc(8)
2686 common /LQGout/ DXSEC,pvalence,q_o,q_s,genproc,genprtyp,
2687 >sch,uch,gproc
2688
2689
2690 double precision CONV
2691 DATA CONV/1.D-9/
2692
2693 DSDXDY=0.
2694 X=xx(1)
2695 Y=xx(2)
2696 Q2=x*y*s
2697 dsdxdy=0.d0
2698 if(Q2.gt.Q2min) then
2699 call LQGDDXS
2700 DSDXDY=(DXSEC(2)+DXSEC(3))*conv
2701 endif
2702 return
2703 end
2704
2705 Subroutine LQGXHMA(XSEC,XSECE,IERR)
2706
2707
2708
2709
2710
2711
2712
2713
2714
2715
2716 implicit none
2717
2718 double precision XSEC,XSECE
2719 integer IERR
2720 external qdens
2721
2722
2723 integer minpts,maxpts,iwk
2724 parameter (minpts=10000)
2725 parameter (maxpts=5000000)
2726 parameter (iwk=1000000)
2727 double precision wk(iwk),a(2),b(2)
2728 double precision eps
2729 parameter (eps=1d-3)
2730 double precision result,relerr
2731 integer nfnevl,ifail
2732 double precision pi
2733 parameter (pi=3.141592653589793d0)
2734
2735
2736 double precision xmax,xmin,ymax,ymin,zmax,zmin,Q2min
2737 common /LQGKinC/ xmax,xmin,ymax,ymin,zmax,zmin,Q2min
2738
2739
2740 double precision Mlq,G1,G2
2741 Integer LQtype,l_o,q_i,q_j
2742 common /LQGproc/ Mlq,G1,G2,LQtype,l_o,q_i,q_j
2743
2744
2745 double precision S,Srad,x,y,z,Q2
2746 common /LQGKinV/ S,Srad,x,y,z,Q2
2747
2748
2749 IERR=0
2750 XSEC=0.d0
2751 XSECE=0.d0
2752 a(1)=xmin
2753 b(1)=xmax
2754 a(2)=ymin
2755 b(2)=ymax
2756
2757 call dadmul(qdens,2,a,b,minpts,maxpts,eps,
2758 + wk,iwk,result,relerr,nfnevl,ifail)
2759
2760 if(ifail.ne.0)then
2761 Write(6,*)'LQGXHMA: Error in Cross Section integration'
2762 IERR=ifail
2763 return
2764 endif
2765 XSEC=S/32./pi*G1*G1*G2*G2/MLQ/MLQ/MLQ/MLQ*result*0.38938
2766 XSECE=relerr*XSEC
2767 return
2768 end
2769
2770 double precision function qdens(n,xx)
2771
2772 implicit none
2773
2774 integer n
2775 double precision xx(*)
2776 double precision UPVs,DNVs,USEAs,DSEAs,STRs,
2777 +CHMs,BOTs,TOPs,GLs
2778 double precision UPVu,DNVu,USEAu,DSEAu,STRu,
2779 +CHMu,BOTu,TOPu,GLu
2780
2781 double precision sh,u
2782 double precision UPQs(3),DNQs(3),UPQBs(3),DNQBs(3)
2783 double precision UPQu(3),DNQu(3),UPQBu(3),DNQBu(3)
2784
2785
2786 double precision xmax,xmin,ymax,ymin,zmax,zmin,Q2min
2787 common /LQGKinC/ xmax,xmin,ymax,ymin,zmax,zmin,Q2min
2788
2789
2790 double precision Mlq,G1,G2
2791 Integer LQtype,l_o,q_i,q_j
2792 common /LQGproc/ Mlq,G1,G2,LQtype,l_o,q_i,q_j
2793
2794
2795 double precision S,Srad,x,y,z,Q2
2796 common /LQGKinV/ S,Srad,x,y,z,Q2
2797
2798
2799 integer echar
2800 double precision ebeam,pbeam
2801 common /LQGbeam/ebeam,pbeam,echar
2802
2803 X=xx(1)
2804 Y=xx(2)
2805
2806 qdens=0.d0
2807 Q2=S*x*y
2808 sh=sqrt(S*X)
2809 u=sqrt(S*X*(1-Y))
2810 if(Q2.gt.Q2min)Then
2811 CALL STRUCTM(X,sh,UPVs
2812 & ,DNVs,USEAs,DSEAs,STRs,CHMs,BOTs,TOPs,GLs)
2813 else
2814 UPVs=0.
2815 DNVs=0.
2816 USEAs=0.
2817 DSEAs=0.
2818 STRs=0.
2819 CHMs=0.
2820 BOTs=0.
2821 TOPs=0.
2822 GLs=0.
2823 endif
2824 if(Q2.gt.Q2min)Then
2825 CALL STRUCTM(X,u,UPVu
2826 & ,DNVu,USEAu,DSEAu,STRu,CHMu,BOTu,TOPu,GLu)
2827 else
2828 UPVu=0.
2829 DNVu=0.
2830 USEAu=0.
2831 DSEAu=0.
2832 STRu=0.
2833 CHMu=0.
2834 BOTu=0.
2835 TOPu=0.
2836 GLu=0.
2837 endif
2838
2839 if(echar.eq.1)then
2840
2841 UPQu(1)=UPVu+USEAu
2842 UPQBu(1)=USEAu
2843 UPQu(2)=CHMu
2844 UPQBu(2)=CHMu
2845 UPQu(3)=TOPu
2846 UPQBu(3)=TOPu
2847 DNQu(1)=DNVu+DSEAu
2848 DNQBu(1)=DSEAu
2849 DNQu(2)=STRu
2850 DNQBu(2)=STRu
2851 DNQu(3)=BOTu
2852 DNQBu(3)=BOTu
2853 elseif(echar.eq.-1)then
2854
2855 UPQu(1)=USEAu
2856 UPQBu(1)=UPVu+USEAu
2857 UPQu(2)=CHMu
2858 UPQBu(2)=CHMu
2859 UPQu(3)=TOPu
2860 UPQBu(3)=TOPu
2861 DNQu(1)=DSEAu
2862 DNQBu(1)=DNVu+DSEAu
2863 DNQu(2)=STRu
2864 DNQBu(2)=STRu
2865 DNQu(3)=BOTu
2866 DNQBu(3)=BOTu
2867 endif
2868
2869 if(echar.eq.1)then
2870
2871 UPQs(1)=UPVs+USEAs
2872 UPQBs(1)=USEAs
2873 UPQs(2)=CHMs
2874 UPQBs(2)=CHMs
2875 UPQs(3)=TOPs
2876 UPQBs(3)=TOPs
2877 DNQs(1)=DNVs+DSEAs
2878 DNQBs(1)=DSEAs
2879 DNQs(2)=STRs
2880 DNQBs(2)=STRs
2881 DNQs(3)=BOTs
2882 DNQBs(3)=BOTs
2883 elseif(echar.eq.-1)then
2884
2885 UPQs(1)=USEAs
2886 UPQBs(1)=UPVs+USEAs
2887 UPQs(2)=CHMs
2888 UPQBs(2)=CHMs
2889 UPQs(3)=TOPs
2890 UPQBs(3)=TOPs
2891 DNQs(1)=DSEAs
2892 DNQBs(1)=DNVs+DSEAs
2893 DNQs(2)=STRs
2894 DNQBs(2)=STRs
2895 DNQs(3)=BOTs
2896 DNQBs(3)=BOTs
2897 endif
2898
2899 if(LQTYPE.eq.1)Then
2900 qdens=UPQBs(q_i)/2.+UPQu(q_j)*(1-y)*(1-y)/2.
2901 elseif(LQTYPE.eq.2)Then
2902 qdens=UPQBs(q_i)/2.+UPQu(q_j)*(1-y)*(1-y)/2.
2903 elseif(LQTYPE.eq.3)Then
2904 qdens=DNQBs(q_i)/2.+DNQu(q_j)*(1-y)*(1-y)/2.
2905 elseif(LQTYPE.eq.4)Then
2906 if(q_j.eq.3)UPQBs(q_i)=0.
2907 if(q_i.eq.3)UPQu(q_j)=0.
2908 qdens=(UPQBs(q_i)+4.*DNQBs(q_i))/2.
2909 > +(UPQu(q_j)+4.*DNQu(q_j))*(1.-Y)*(1.-Y)/2.
2910 elseif(LQTYPE.eq.5)Then
2911 qdens=DNQBs(q_i)*(1-y)*(1-y)*2.+DNQu(q_j)*2.
2912 elseif(LQTYPE.eq.6)Then
2913 if(q_j.eq.3)UPQBs(q_i)=0.
2914 if(q_i.eq.3)UPQu(q_j)=0.
2915 qdens=(UPQBs(q_i)+DNQBs(q_i))*(1-y)*(1-y)*2.
2916 > +(UPQu(q_j)+DNQu(q_j))*2.
2917 elseif(LQTYPE.eq.7)Then
2918 qdens=UPQBs(q_i)*(1-y)*(1-y)*2.+UPQu(q_j)*2.
2919 elseif(LQTYPE.eq.8)Then
2920 qdens=DNQs(q_i)*2.*(1.-Y)*(1.-Y)+DNQBu(q_j)*2.
2921 elseif(LQTYPE.eq.9)Then
2922 qdens=DNQs(q_i)*2.*(1.-Y)*(1.-Y)+DNQBu(q_j)*2.
2923 elseif(LQTYPE.eq.10)Then
2924 qdens=UPQs(q_i)*2.*(1.-Y)*(1.-Y)+UPQBu(q_j)*2.
2925 elseif(LQTYPE.eq.11)Then
2926 if(q_j.eq.3)UPQs(q_i)=0.
2927 if(q_i.eq.3)UPQBu(q_j)=0.
2928 qdens=(4.*UPQs(q_i)+DNQs(q_i))*2.*(1.-Y)*(1.-Y)
2929 > +(4.*UPQBu(q_j)+DNQBu(q_j))*2.
2930 elseif(LQTYPE.EQ.12)Then
2931 qdens=UPQs(q_i)/2.+UPQBu(q_j)*(1.-Y)*(1.-Y)/2.
2932 elseif(LQTYPE.EQ.13)Then
2933 if(q_j.eq.3)UPQs(q_i)=0.
2934 if(q_i.eq.3)UPQBu(q_j)=0.
2935 qdens=(UPQs(q_i)+DNQs(q_i))/2.
2936 > +(UPQBu(q_j)+DNQBu(q_j))*(1.-Y)*(1.-Y)/2.
2937 elseif(LQTYPE.EQ.14)Then
2938 qdens=DNQs(q_i)/2.+DNQBu(q_j)*(1.-Y)*(1.-Y)/2.
2939 endif
2940 return
2941 end
2942