updating apd and fee parameters

[u/mrichter/AliRoot.git] / PYTHIA6 / AliPythia.cxx

diff --git a/PYTHIA6/AliPythia.cxx b/PYTHIA6/AliPythia.cxx
index 6520a1350466c364ae90561aac3f688493b89978..b31c9481445643248a9081dee20a176022a13a2a 100644 (file)
--- a/PYTHIA6/AliPythia.cxx
+++ b/PYTHIA6/AliPythia.cxx
@@ -18,9 +18,10 @@
 
 #include "AliPythia.h"
 #include "AliPythiaRndm.h"
-#include "../FASTSIM/AliFastGlauber.h"
-#include "../FASTSIM/AliQuenchingWeights.h"
+#include "AliFastGlauber.h"
+#include "AliQuenchingWeights.h"
 #include "TVector3.h"
+#include "PyquenCommon.h"
 
 ClassImp(AliPythia)
 
@@ -58,6 +59,8 @@ AliPythia::AliPythia():
     fStrucFunc(kCTEQ5L),
     fXJet(0.),
     fYJet(0.),
+    fNGmax(30),
+    fZmax(0.97),
     fGlauber(0),
     fQuenchingWeights(0)
 {
@@ -78,6 +81,8 @@ AliPythia::AliPythia(const AliPythia& pythia):
     fStrucFunc(kCTEQ5L),
     fXJet(0.),
     fYJet(0.),
+    fNGmax(30),
+    fZmax(0.97),
     fGlauber(0),
     fQuenchingWeights(0)
 {
@@ -95,18 +100,18 @@ void AliPythia::ProcInit(Process_t process, Float_t energy, StrucFunc_t strucfun
     fEcms = energy;
     fStrucFunc = strucfunc;
 //...Switch off decay of pi0, K0S, Lambda, Sigma+-, Xi0-, Omega-.
-    SetMDCY(Pycomp(111) ,1,0);
-    SetMDCY(Pycomp(310) ,1,0);
-    SetMDCY(Pycomp(3122),1,0);
-    SetMDCY(Pycomp(3112),1,0);
-    SetMDCY(Pycomp(3212),1,0);
-    SetMDCY(Pycomp(3222),1,0);
-    SetMDCY(Pycomp(3312),1,0);
-    SetMDCY(Pycomp(3322),1,0);
-    SetMDCY(Pycomp(3334),1,0);
+    SetMDCY(Pycomp(111) ,1,0); // pi0
+    SetMDCY(Pycomp(310) ,1,0); // K0S
+    SetMDCY(Pycomp(3122),1,0); // kLambda
+    SetMDCY(Pycomp(3112),1,0); // sigma -
+    SetMDCY(Pycomp(3212),1,0); // sigma 0 
+    SetMDCY(Pycomp(3222),1,0); // sigma +
+    SetMDCY(Pycomp(3312),1,0); // xi - 
+    SetMDCY(Pycomp(3322),1,0); // xi 0
+    SetMDCY(Pycomp(3334),1,0); // omega-
     // Select structure function 
     SetMSTP(52,2);
-    SetMSTP(51,strucfunc);
+    SetMSTP(51, AliStructFuncType::PDFsetIndex(strucfunc));
     // Particles produced in string fragmentation point directly to either of the two endpoints
     // of the string (depending in the side they were generated from).
     SetMSTU(16,2);
@@ -260,7 +265,7 @@ void AliPythia::ProcInit(Process_t process, Float_t energy, StrucFunc_t strucfun
        SetMSUB(94,1);             // double diffraction
        SetMSUB(95,1);             // low pt production
 
-        SetMSTP(51,kCTEQ6ll);      // CTEQ6ll pdf
+        SetMSTP(51,AliStructFuncType::PDFsetIndex(kCTEQ6ll));      // CTEQ6ll pdf
        SetMSTP(52,2);
        SetMSTP(68,1);
        SetMSTP(70,2);
@@ -305,7 +310,7 @@ void AliPythia::ProcInit(Process_t process, Float_t energy, StrucFunc_t strucfun
        SetMSEL(1);
  // Pythia Tune A (CDF)
  //
-       SetPARP(67,4.);            // Regulates Initial State Radiation
+       SetPARP(67,2.5);           // Regulates Initial State Radiation (value from best fit to D0 dijet analysis)
        SetMSTP(82,4);             // Double Gaussian Model
        SetPARP(82,2.0);           // [GeV]    PT_min at Ref. energy
        SetPARP(84,0.4);           // Core radius
@@ -546,9 +551,8 @@ void AliPythia::ProcInit(Process_t process, Float_t energy, StrucFunc_t strucfun
 //
 //  Initialize PYTHIA
     SetMSTP(41,1);   // all resonance decays switched on
-
     Initialize("CMS","p","p",fEcms);
-
+    
 }
 
 Int_t AliPythia::CheckedLuComp(Int_t kf)
@@ -657,7 +661,7 @@ void AliPythia::Pyrobo(Int_t imi, Int_t ima, Double_t the, Double_t phi, Double_
 
 
 
-void AliPythia::InitQuenching(Float_t cMin, Float_t cMax, Float_t k, Int_t iECMethod)
+void AliPythia::InitQuenching(Float_t cMin, Float_t cMax, Float_t k, Int_t iECMethod, Float_t zmax, Int_t ngmax)
 {
 // Initializes 
 // (1) The quenching model using quenching weights according to C. Salgado and U. Wiedemann
@@ -672,6 +676,9 @@ void AliPythia::InitQuenching(Float_t cMin, Float_t cMax, Float_t k, Int_t iECMe
     fQuenchingWeights->InitMult();
     fQuenchingWeights->SetK(k);
     fQuenchingWeights->SetECMethod(AliQuenchingWeights::kECMethod(iECMethod));
+    fNGmax = ngmax;
+    fZmax  = zmax;
+    
 }
 
 
@@ -789,7 +796,7 @@ void  AliPythia::Quench()
            //
            // Avoid complete loss
            //
-           if (fZQuench[j] == 1.) fZQuench[j] = 0.97;
+           if (fZQuench[j] == 1.) fZQuench[j] = fZmax;
            //
            // Some debug printing
 
@@ -818,7 +825,7 @@ void  AliPythia::Quench()
        if (!quenched[isys]) continue;
        
        nGluon[isys]   = 1 + Int_t(fZQuench[isys] / (1. - fZQuench[isys]));
-       if (nGluon[isys] > 30) nGluon[isys] = 30;
+       if (nGluon[isys] > fNGmax) nGluon[isys] = fNGmax;
        zquench[isys] = 1. - TMath::Power(1. - fZQuench[isys], 1./Double_t(nGluon[isys]));
        wjtKick[isys]  = wjtKick[isys] / TMath::Sqrt(Double_t(nGluon[isys]));
 
@@ -954,9 +961,13 @@ void  AliPythia::Quench()
                    //
                    //     Calculate new px, py
                    //
-                   Double_t pxNew   = jtNew / jt * pxs;
-                   Double_t pyNew   = jtNew / jt * pys;        
+                   Double_t pxNew   = 0;
+                   Double_t pyNew   = 0;
                    
+                   if (jt>0) {
+                     pxNew = jtNew / jt * pxs;
+                     pyNew = jtNew / jt * pys;
+                   }   
 //                 Double_t dpx = pxs - pxNew;
 //                 Double_t dpy = pys - pyNew;
 //                 Double_t dpz = pl  - plNew;
@@ -1203,9 +1214,25 @@ void  AliPythia::Quench()
 void AliPythia::Pyquen(Double_t a, Int_t ibf, Double_t b)
 {
     // Igor Lokthine's quenching routine
+    // http://lokhtin.web.cern.ch/lokhtin/pyquen/pyquen.txt
+
     pyquen(a, ibf, b);
 }
 
+void AliPythia::SetPyquenParameters(Double_t t0, Double_t tau0, Int_t nf, Int_t iengl, Int_t iangl)
+{
+    // Set the parameters for the PYQUEN package.
+    // See comments in PyquenCommon.h
+    
+    
+    PYQPAR.t0    = t0;
+    PYQPAR.tau0  = tau0;
+    PYQPAR.nf    = nf;
+    PYQPAR.iengl = iengl;
+    PYQPAR.iangl = iangl;
+}
+
+
 void AliPythia::Pyevnw()
 {
     // New multiple interaction scenario
@@ -1250,7 +1277,7 @@ void AliPythia::AtlasTuning()
 {
     //
     // Configuration for the ATLAS tuning
-        SetMSTP(51, kCTEQ5L);      // CTEQ5L pdf
+        SetMSTP(51, AliStructFuncType::PDFsetIndex(kCTEQ5L));      // CTEQ5L pdf
        SetMSTP(81,1);             // Multiple Interactions ON
        SetMSTP(82,4);             // Double Gaussian Model
        SetPARP(81,1.9);           // Min. pt for multiple interactions (default in 6.2-14)