C...Double precision and integer declarations.
IMPLICIT DOUBLE PRECISION(A-H, O-Z)
IMPLICIT INTEGER(I-N)
- INTEGER PYK,PYCHGE,PYCOMP
+C INTEGER PYK,PYCHGE,PYCOMP
C...Commonblocks.
COMMON/PYDAT1/MSTU(200),PARU(200),MSTJ(200),PARJ(200)
COMMON/PYDAT2/KCHG(500,4),PMAS(500,4),PARF(2000),VCKM(4,4)
CHARACTER CHFRAM*12,CHBEAM*12,CHTARG*12,CHLH(2)*6
C...Interface to PDFLIB.
- COMMON/W50512/QCDL4,QCDL5
- SAVE /W50512/
+ COMMON/LW50512/QCDL4,QCDL5
+ SAVE /LW50512/
DOUBLE PRECISION VALUE(20),QCDL4,QCDL5
CHARACTER*20 PARM(20)
DATA VALUE/20*0D0/,PARM/20*' '/
WRITE(MSTU(11),5100)
STOP
ENDIF
-
C...Initial values for some counters.
N=0
MINT(5)=MINT(5)+1
MINT(51)=0
C...Evaluate suppression factors due to non-simulated channels.
- IF(KCHG(KC,3).EQ.0) THEN
- WIDS(KC,1)=((WDTE(0,1)+WDTE(0,2))**2+
- & 2D0*(WDTE(0,1)+WDTE(0,2))*(WDTE(0,4)+WDTE(0,5))+
- & 2D0*WDTE(0,4)*WDTE(0,5))/WDTP(0)**2
- WIDS(KC,2)=(WDTE(0,1)+WDTE(0,2)+WDTE(0,4))/WDTP(0)
- WIDS(KC,3)=0D0
- WIDS(KC,4)=0D0
- WIDS(KC,5)=0D0
- ELSE
- IF(MWID(KC).EQ.3) MINT(63)=1
- CALL PYWIDT(-KF,PMR**2,WDTPM,WDTEM)
- MINT(51)=0
- WIDS(KC,1)=((WDTE(0,1)+WDTE(0,2))*(WDTEM(0,1)+WDTEM(0,3))+
- & (WDTE(0,1)+WDTE(0,2))*(WDTEM(0,4)+WDTEM(0,5))+
- & (WDTE(0,4)+WDTE(0,5))*(WDTEM(0,1)+WDTEM(0,3))+
- & WDTE(0,4)*WDTEM(0,5)+WDTE(0,5)*WDTEM(0,4))/WDTP(0)**2
- WIDS(KC,2)=(WDTE(0,1)+WDTE(0,2)+WDTE(0,4))/WDTP(0)
- WIDS(KC,3)=(WDTEM(0,1)+WDTEM(0,3)+WDTEM(0,4))/WDTP(0)
- WIDS(KC,4)=((WDTE(0,1)+WDTE(0,2))**2+
- & 2D0*(WDTE(0,1)+WDTE(0,2))*(WDTE(0,4)+WDTE(0,5))+
- & 2D0*WDTE(0,4)*WDTE(0,5))/WDTP(0)**2
- WIDS(KC,5)=((WDTEM(0,1)+WDTEM(0,3))**2+
- & 2D0*(WDTEM(0,1)+WDTEM(0,3))*(WDTEM(0,4)+WDTEM(0,5))+
- & 2D0*WDTEM(0,4)*WDTEM(0,5))/WDTP(0)**2
+C...AM
+C...Protection against division by 0 since rho_21_tc is causing problem here
+ IF (WDTP(0) .GT. 0.) THEN
+
+ IF(KCHG(KC,3).EQ.0) THEN
+ WIDS(KC,1)=((WDTE(0,1)+WDTE(0,2))**2+
+ & 2D0*(WDTE(0,1)+WDTE(0,2))*(WDTE(0,4)+WDTE(0,5))+
+ & 2D0*WDTE(0,4)*WDTE(0,5))/WDTP(0)**2
+ WIDS(KC,2)=(WDTE(0,1)+WDTE(0,2)+WDTE(0,4))/WDTP(0)
+ WIDS(KC,3)=0D0
+ WIDS(KC,4)=0D0
+ WIDS(KC,5)=0D0
+ ELSE
+ IF(MWID(KC).EQ.3) MINT(63)=1
+ CALL PYWIDT(-KF,PMR**2,WDTPM,WDTEM)
+ MINT(51)=0
+ WIDS(KC,1)=((WDTE(0,1)+WDTE(0,2))*(WDTEM(0,1)+WDTEM(0,3))+
+ & (WDTE(0,1)+WDTE(0,2))*(WDTEM(0,4)+WDTEM(0,5))+
+ & (WDTE(0,4)+WDTE(0,5))*(WDTEM(0,1)+WDTEM(0,3))+
+ & WDTE(0,4)*WDTEM(0,5)+WDTE(0,5)*WDTEM(0,4))/WDTP(0)**2
+ WIDS(KC,2)=(WDTE(0,1)+WDTE(0,2)+WDTE(0,4))/WDTP(0)
+ WIDS(KC,3)=(WDTEM(0,1)+WDTEM(0,3)+WDTEM(0,4))/WDTP(0)
+ WIDS(KC,4)=((WDTE(0,1)+WDTE(0,2))**2+
+ & 2D0*(WDTE(0,1)+WDTE(0,2))*(WDTE(0,4)+WDTE(0,5))+
+ & 2D0*WDTE(0,4)*WDTE(0,5))/WDTP(0)**2
+ WIDS(KC,5)=((WDTEM(0,1)+WDTEM(0,3))**2+
+ & 2D0*(WDTEM(0,1)+WDTEM(0,3))*(WDTEM(0,4)+WDTEM(0,5))+
+ & 2D0*WDTEM(0,4)*WDTEM(0,5))/WDTP(0)**2
+ ENDIF
+
ENDIF
-
C...Set resonance widths and branching ratios;
C...also on/off switch for decays.
IF(MWID(KC).EQ.1.OR.MWID(KC).EQ.3) THEN
&XPPI(-6:6),XPPR(-6:6)
C...Interface to PDFLIB.
- COMMON/W50513/XMIN,XMAX,Q2MIN,Q2MAX
- SAVE /W50513/
+ COMMON/LW50513/XMIN,XMAX,Q2MIN,Q2MAX
+ SAVE /LW50513/
DOUBLE PRECISION XX,QQ,UPV,DNV,USEA,DSEA,STR,CHM,BOT,TOP,GLU,
&VALUE(20),XMIN,XMAX,Q2MIN,Q2MAX
CHARACTER*20 PARM(20)
DIMENSION XPEL(-25:25),XPGA(-6:6),SXP(0:6)
C...Interface to PDFLIB.
- COMMON/W50513/XMIN,XMAX,Q2MIN,Q2MAX
- SAVE /W50513/
+ COMMON/LW50513/XMIN,XMAX,Q2MIN,Q2MAX
+ SAVE /LW50513/
DOUBLE PRECISION XX,QQ,UPV,DNV,USEA,DSEA,STR,CHM,BOT,TOP,GLU,
&VALUE(20),XMIN,XMAX,Q2MIN,Q2MAX
CHARACTER*20 PARM(20)
SAVE /PYJETS/,/PYDAT1/,/PYDAT2/,/PYDAT3/,/PYINT4/
C...Local array.
DIMENSION PS(2,6),IJOIN(100)
-
C...Initialize and reset.
MSTU(24)=0
IF(MSTU(12).GE.1) CALL PYLIST(0)
KF=ISIGN(1000*KFLD+100*KFLE+10*KFLF+KFLS,KFL1)
IF(KSIG.EQ.0) KF=ISIGN(1000*KFLD+100*KFLF+10*KFLE+KFLS,KFL1)
ENDIF
+C -------------------------------------------------------------------------
+C Extracted from a private e-mail exchange with Torbjorn Sjostrand
+C
+C No, Lambda(1520) is not included and not foreseen.
+C So if you want it in Pythia, it would have to be a hack.
+C What you could do is:
+C 1) In PYKFDI, just before the RETURN above label 140, you could check if
+C a Lambda, Sigma0 or Sigma*0 has been produced, and with some small
+C probability switch such a particle to the Lambda(1520) code. That is,
+C if KF = 3122, 3212, or 3214 and a random number below some number, switch
+C to KF = 3124. (And correspondingly for anticparticles.)
+C 2) Use the PYUPDA routine (see manual) to include particle and decay data
+C for the Lambda(1520).
+C -------------------------------------------------------------------------
+
+C IF (IABS(KF).EQ.3122) THEN
+C Converting a fraction (0.20) of Lambda0 to Lambda(1520) + c.c.
+C This fraction is based on the experimental measurement at ISR
+C Bobbink 83, NP B217,11 (1983)
+C The region 0.5 < XF < 1.0 has been extrapolated to XF=0
+C IF(PYR(0).LE.0.20) KF=ISIGN(3124,KF)
+C ENDIF
+
+C IF(IABS(KF).EQ.3212) THEN
+C Converting a fraction (0.20) of Sigma0 to Lambda(1520) + c.c.
+C We suppose the same fraction as for Lambda0
+C IF(PYR(0).LE.0.20) KF=ISIGN(3124,KF)
+C ENDIF
+
+C IF (IABS(KF).EQ.3214) THEN
+C Converting a fraction (0.30) of Sigma0(1385) to Lambda(1520) + c.c.
+C This is conservative extimate supposing that the ratio
+C scales as (M_Sigma1385/M_Lambda0)^2 ~ 1.5
+C IF(PYR(0).LE.0.30) KF=ISIGN(3124,KF)
+C ENDIF
RETURN
C...Use tabulated probabilities to select new flavour and hadron.
C...Integral of Altarelli-Parisi z kernel for QCD.
C...(Includes squark and gluino; with factor N_C/C_F extra for latter).
+
+ FMED = PARJ(200)
IF(MSTJ(49).EQ.0.AND.KFL(1).EQ.21) THEN
- FBR=6D0*LOG((1D0-ZC)/ZC)+MSTJ(45)*0.5D0
+C Nestor
+ FBR=(1.D0+FMED)*6D0*LOG((1D0-ZC)/ZC)+MSTJ(45)*0.5D0
ELSEIF(MSTJ(49).EQ.0) THEN
- FBR=(8D0/3D0)*LOG((1D0-ZC)/ZC)
+C Nestor
+ FBR=(1.D0+FMED)*(8D0/3D0)*LOG((1D0-ZC)/ZC)
IF(IGLUI.EQ.1.AND.IR.GE.31) FBR=FBR*(9D0/4D0)
C...Integral of Altarelli-Parisi z kernel for scalar gluon.
C...Integral of Altarelli-Parisi z kernel for Abelian vector gluon.
ELSEIF(KFL(1).EQ.21) THEN
- FBR=6D0*MSTJ(45)*(0.5D0-ZC)
+ FBR=(1.D0+FMED)*6D0*MSTJ(45)*(0.5D0-ZC)
ELSE
- FBR=2D0*LOG((1D0-ZC)/ZC)
+ FBR=(1.D0+FMED)*2D0*LOG((1D0-ZC)/ZC)
ENDIF
C...Reset QCD probability for colourless.
150 CONTINUE
ENDIF
C...Remove duplicate old decay data.
- IF(KCREP.NE.0.AND.MDCY(KCREP,3).GT.0) THEN
- IDCREP=MDCY(KCREP,2)
- NDCREP=MDCY(KCREP,3)
- DO 160 I=1,KCC
- IF(MDCY(I,2).GT.IDCREP) MDCY(I,2)=MDCY(I,2)-NDCREP
- 160 CONTINUE
- DO 180 I=IDCREP,NDC-NDCREP
- MDME(I,1)=MDME(I+NDCREP,1)
- MDME(I,2)=MDME(I+NDCREP,2)
- BRAT(I)=BRAT(I+NDCREP)
- DO 170 J=1,5
- KFDP(I,J)=KFDP(I+NDCREP,J)
- 170 CONTINUE
- 180 CONTINUE
- NDC=NDC-NDCREP
- KC=KCREP
- ELSEIF(KCREP.NE.0) THEN
- KC=KCREP
+ IF(KCREP.NE.0) THEN
+ IF(MDCY(KCREP,3).GT.0) THEN
+ IDCREP=MDCY(KCREP,2)
+ NDCREP=MDCY(KCREP,3)
+ DO 160 I=1,KCC
+ IF(MDCY(I,2).GT.IDCREP) MDCY(I,2)=MDCY(I,2)-NDCREP
+ 160 CONTINUE
+ DO 180 I=IDCREP,NDC-NDCREP
+ MDME(I,1)=MDME(I+NDCREP,1)
+ MDME(I,2)=MDME(I+NDCREP,2)
+ BRAT(I)=BRAT(I+NDCREP)
+ DO 170 J=1,5
+ KFDP(I,J)=KFDP(I+NDCREP,J)
+ 170 CONTINUE
+ 180 CONTINUE
+ NDC=NDC-NDCREP
+ KC=KCREP
+ ELSE
+ KC=KCREP
+ ENDIF
ELSE
KCC=KCC+1
KC=KCC
C IDATI(6)=ISEC
C...Example 4: GNU LINUX libU77, SunOS.
- CALL IDATE(IDTEMP)
- IDATI(1)=IDTEMP(3)
- IDATI(2)=IDTEMP(2)
- IDATI(3)=IDTEMP(1)
- CALL ITIME(IDTEMP)
- IDATI(4)=IDTEMP(1)
- IDATI(5)=IDTEMP(2)
- IDATI(6)=IDTEMP(3)
+c CALL IDATE(IDTEMP)
+c IDATI(1)=IDTEMP(3)
+c IDATI(2)=IDTEMP(2)
+c IDATI(3)=IDTEMP(1)
+c CALL ITIME(IDTEMP)
+c IDATI(4)=IDTEMP(1)
+c IDATI(5)=IDTEMP(2)
+c IDATI(6)=IDTEMP(3)
C...Common code to ensure right century.
IDATI(1)=2000+MOD(IDATI(1),100)
git://git.uio.no / u / mrichter / AliRoot.git / blobdiff