Analysis Software
Documentation for sPHENIX simulation software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
luxkfl.f
Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file luxkfl.f
1 
2 C*********************************************************************
3 
4  SUBROUTINE luxkfl(KFL,ECM,ECMC,KFLC)
5 
6 C...Purpose: to select flavour for produced qqbar pair.
7  common/ludat1/mstu(200),paru(200),mstj(200),parj(200)
8  SAVE /ludat1/
9  common/ludat2/kchg(500,3),pmas(500,4),parf(2000),vckm(4,4)
10  SAVE /ludat2/
11 
12 C...Calculate maximum weight in QED or QFD case.
13  IF(mstj(102).LE.1) THEN
14  rfmax=4./9.
15  ELSE
16  poll=1.-parj(131)*parj(132)
17  sff=1./(16.*paru(102)*(1.-paru(102)))
18  sfw=ecmc**4/((ecmc**2-parj(123)**2)**2+(parj(123)*parj(124))**2)
19  sfi=sfw*(1.-(parj(123)/ecmc)**2)
20  ve=4.*paru(102)-1.
21  hf1i=sfi*sff*(ve*poll+parj(132)-parj(131))
22  hf1w=sfw*sff**2*((ve**2+1.)*poll+2.*ve*(parj(132)-parj(131)))
23  rfmax=max(4./9.*poll-4./3.*(1.-8.*paru(102)/3.)*hf1i+
24  & ((1.-8.*paru(102)/3.)**2+1.)*hf1w,1./9.*poll+2./3.*
25  & (-1.+4.*paru(102)/3.)*hf1i+((-1.+4.*paru(102)/3.)**2+1.)*hf1w)
26  ENDIF
27 
28 C...Choose flavour. Gives charge and velocity.
29  ntry=0
30  100 ntry=ntry+1
31  IF(ntry.GT.100) THEN
32  CALL luerrm(14,'(LUXKFL:) caught in an infinite loop')
33  kflc=0
34  RETURN
35  ENDIF
36  kflc=kfl
37  IF(kfl.LE.0) kflc=1+int(mstj(104)*rlu(0))
38  mstj(93)=1
39  pmq=ulmass(kflc)
40  IF(ecm.LT.2.*pmq+parj(127)) goto 100
41  qf=kchg(kflc,1)/3.
42  vq=1.
43  IF(mod(mstj(103),2).EQ.1) vq=sqrt(max(0.,1.-(2.*pmq/ecmc)**2))
44 
45 C...Calculate weight in QED or QFD case.
46  IF(mstj(102).LE.1) THEN
47  rf=qf**2
48  rfv=0.5*vq*(3.-vq**2)*qf**2
49  ELSE
50  vf=sign(1.,qf)-4.*qf*paru(102)
51  rf=qf**2*poll-2.*qf*vf*hf1i+(vf**2+1.)*hf1w
52  rfv=0.5*vq*(3.-vq**2)*(qf**2*poll-2.*qf*vf*hf1i+vf**2*hf1w)+
53  & vq**3*hf1w
54  ENDIF
55 
56 C...Weighting or new event (radiative photon). Cross-section update.
57  IF(kfl.LE.0.AND.rf.LT.rlu(0)*rfmax) goto 100
58  parj(158)=parj(158)+1.
59  IF(ecmc.LT.2.*pmq+parj(127).OR.rfv.LT.rlu(0)*rf) kflc=0
60  IF(mstj(107).LE.0.AND.kflc.EQ.0) goto 100
61  IF(kflc.NE.0) parj(159)=parj(159)+1.
62  parj(144)=parj(157)*parj(159)/parj(158)
63  parj(148)=parj(144)*86.8/ecm**2
64 
65  RETURN
66  END