add histograms with sum of split energy or pT of split clusters after selection as pi0

[u/mrichter/AliRoot.git] / PYTHIA6 / pythia6.4.25 / pythia-6.4.25.f

diff --git a/PYTHIA6/pythia6.4.25/pythia-6.4.25.f b/PYTHIA6/pythia6.4.25/pythia-6.4.25.f
index b888a24d3c8907800ef6e8a720f2bb9bc481ad07..cbaef1e55ca85624c9a5b9b17afab03c606c2f2e 100644 (file)
--- a/PYTHIA6/pythia6.4.25/pythia-6.4.25.f
+++ b/PYTHIA6/pythia6.4.25/pythia-6.4.25.f
@@ -3583,6 +3583,11 @@ C...Stop if no subprocesses on.
         CALL PYSTOP(1)
       ENDIF
  
+      DO 2 I = 1, 4000
+         DO 1 J = 1, 5
+            V(I,J) = 0.
+ 1       ENDDO
+ 2    ENDDO
 C...Initial values for some counters.
       MSTU(1)=0
       MSTU(2)=0
@@ -15309,8 +15314,12 @@ C...Evolve q -> q gamma separately, pick it if larger pT.
 C...  Evolve massive quark creation separately.
         MCRQQ=0
         IF (MQMASS.NE.0) THEN
-          PT2CR=(RMQ2+VINT(18))*(RML**(TPM/(TPL*PYR(0)**(-TML/WN)-TPM)))
-     &         -VINT(18)
+           IF (WN .EQ. 0.) THEN
+              ARG = -1.
+           ELSE
+              ARG = TPM/(TPL*PYR(0)**(-TML/WN)-TPM)
+           ENDIF
+           PT2CR=(RMQ2+VINT(18))*(RML**ARG)-VINT(18)
 C...If massive quark also on opposite side, ensure sufficient remaining 
 C...phase space also for creation of that quark
           TMINQQ = TMIN
@@ -24080,7 +24089,7 @@ C...  Total string piece energy, momentum squared, and components
      &           + (P(ICR(IS),3)+ P(ICR(IS+1),3))**2
             PZTS = P(ICR(IS),1)*TX+P(ICR(IS),2)*TY+P(ICR(IS),3)*TZ 
      &           + P(ICR(IS+1),1)*TX+P(ICR(IS+1),2)*TY+P(ICR(IS+1),3)*TZ
-            PTTS = SQRT(PPS2 - PZTS**2)
+            PTTS = SQRT(MAX(0.D0,PPS2 - PZTS**2))
 C...  Mass of string piece in units of mpi (at least 1)
             RMPI2  = 0.135D0 
             RM2STR = MAX(RMPI2,EES**2 - PPS2)
@@ -69214,6 +69223,43 @@ C.. Form baryon. Distinguish Lambda- and Sigmalike baryons.
       ENDIF
       RETURN
  
+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 -------------------------------------------------------------------------
+ 
+      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
+         IF(PYR(0).LE.0.20) KF=ISIGN(3124,KF)
+      ENDIF
+
+      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
+         IF(PYR(0).LE.0.20) KF=ISIGN(3124,KF)
+      ENDIF
+
+      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 
+         IF(PYR(0).LE.0.30) KF=ISIGN(3124,KF)
+      ENDIF
+      RETURN
+
 C...Use tabulated probabilities to select new flavour and hadron.
   140 IF(KTAB2.EQ.0.AND.MSTJ(12).LE.0) THEN
         KT3L=1
@@ -80404,8 +80450,14 @@ C...
      >        XMIN,XMAX,Q2MIN,Q2MAX)
          IF (MSTP(194) .EQ. 0) THEN 
             CALL SETLHAPARM("EKS98")
-         ELSE
+         ELSE IF (MSTP(194) .EQ.  9) THEN
+            CALL SETLHAPARM("EPS09LO")
+         ELSE IF (MSTP(194) .EQ. 19) THEN
+            CALL SETLHAPARM("EPS09NLO")
+         ELSE IF (MSTP(194) .EQ.  8) THEN
             CALL SETLHAPARM("EPS08")
+         ELSE
+            CALL SETLHAPARM("EPS09LO")
          ENDIF
       ELSE 
          write(6,*) "-> pdfset"