From Theodor Rascanu: These changes are necessary since we need to have the code...
authoragheata <agheata@f7af4fe6-9843-0410-8265-dc069ae4e863>
Tue, 17 Jul 2012 08:44:18 +0000 (08:44 +0000)
committeragheata <agheata@f7af4fe6-9843-0410-8265-dc069ae4e863>
Tue, 17 Jul 2012 08:44:18 +0000 (08:44 +0000)
EVGEN/AliGenEMCocktail.cxx
EVGEN/AliGenEMlib.cxx
EVGEN/AliGenEMlib.h
EVGEN/AliGenLib.h
EVGEN/AliGenParam.cxx
EVGEN/AliGenParam.h
STEER/STEER/AliGenerator.cxx
STEER/STEER/AliGenerator.h

index e1a5efc..2c78959 100644 (file)
@@ -193,6 +193,9 @@ void AliGenEMCocktail::Generate()
   
 // Loop over generators and generate events
   Int_t igen = 0;
+  Float_t evPlane;
+  Rndm(&evPlane,1);
+  evPlane*=TMath::Pi()*2;
   while((entry = (AliGenCocktailEntry*)next())) {
     gen = entry->Generator();
     gen->SetVertex(fVertex.At(0), fVertex.At(1), fVertex.At(2));
@@ -202,6 +205,7 @@ void AliGenEMCocktail::Generate()
       igen++;  
       if (igen == 1) entry->SetFirst(0);               
       else  entry->SetFirst((partArray->GetEntriesFast())+1);
+      gen->SetEventPlane(evPlane);
       gen->SetNumberParticles(fNPart);         
       gen->Generate();
       entry->SetLast(partArray->GetEntriesFast());
index 5ae592f..1851e87 100644 (file)
 
 ClassImp(AliGenEMlib)
 
+AliGenEMlib::Centbins_t AliGenEMlib::fCentbin=k2030;
+
+Double_t AliGenEMlib::CrossOverLc(const double a, const double b, const double x){
+if(x<b-a/2) return 1.0;
+  else if(x>b+a/2) return 0.0;
+  else return cos(((x-b)/a+0.5)*TMath::Pi())/2+0.5;
+}
+Double_t AliGenEMlib::CrossOverRc(const double a, const double b, const double x){
+  return 1-CrossOverLc(a,b,x);
+}
+
+const Double_t AliGenEMlib::fpTparam[8][14] = {
+  // charged pion
+  { 0.0005650000, 4.3540000000, 7.0070000000, 0.2472, 0.0000000000, 1.0000000000, 0.0005650000, 4.3540000000, 7.0070000000, 0.2472, 9999.0000000, 1.0000000000, 0.0000000000, 0.0000000000 },   //
+  { 0.0005650000, 4.3540000000, 7.0070000000, 0.2472, 0.0000000000, 1.0000000000, 0.0005650000, 4.3540000000, 7.0070000000, 0.2472, 9999.0000000, 1.0000000000, 0.0000000000, 0.0000000000 },   //
+  { 3.217746e+03, 1.632570e+00, 9.340162e+00, 1.0000, 2.650000e+00, 4.200635e+00, 9.039589e+08, 1.191548e-01, 6.907293e+00, 1.0000, 6.750000e+00, 4.099970e+00, 6.036772e+01, 5.928279e+00 },   //10-20
+  { 2.125480e+03, 1.711882e+00, 9.585665e+00, 1.0000, 3.100000e+00, 5.041511e+00, 2.431808e+08, 1.155071e-01, 6.574663e+00, 1.0000, 6.250000e+00, 1.842070e+00, 3.928902e+01, 5.860970e+00 },   //20-30
+  { 1.577897e+03, 1.411948e+00, 8.638815e+00, 1.0000, 2.550000e+00, 4.066432e+00, 3.774439e+08, 1.000330e-01, 6.535971e+00, 1.0000, 6.750000e+00, 2.482514e+00, 3.495494e+01, 5.954321e+00 },   //30-40
+  { 7.061859e+02, 1.223810e+00, 8.532113e+00, 1.0000, 1.350000e+00, 1.956213e+00, 1.318169e+04, 5.658401e-01, 7.157575e+00, 1.0000, 5.250000e+00, 3.900000e+00, 1.958841e+01, 6.114787e+00 },   //40-50
+  { 7.061859e+02, 1.223810e+00, 8.532113e+00, 1.0000, 1.350000e+00, 1.956213e+00, 1.318169e+04, 5.658401e-01, 7.157575e+00, 1.0000, 5.250000e+00, 3.900000e+00, 1.958841e+01, 6.114787e+00 },   //
+  { 0.0005650000, 4.3540000000, 7.0070000000, 0.2472, 0.0000000000, 1.0000000000, 0.0005650000, 4.3540000000, 7.0070000000, 0.2472, 9999.0000000, 1.0000000000, 0.0000000000, 0.0000000000 } }; //pp
+
+const Double_t AliGenEMlib::fv2param[2][8][9] = { { 
+    // charged pion
+    {  3.971545e-02, 2.111261e+00, 6.499885e-02, 3.659836e+00, 7.812234e+00, 1.479471e-02, 1.241708e+00, 2.817950e+00, 1.910663e-01 }   //0-5
+    ,{ 7.424082e-02, 2.417588e+00, 7.891084e-02, 2.232841e+00, 3.398147e+00, 3.410206e-02, 4.138276e-01,-1.152430e+00, 5.729093e-01 }   //5-10
+    ,{ 1.094608e-01, 2.357420e+00, 7.614904e-02, 2.294245e+00, 3.538364e+00, 4.932739e-02, 4.557926e-01, 9.218702e-03, 6.428540e-01 }   //10-20
+    ,{ 1.456377e-01, 2.322408e+00, 7.747166e-02, 2.148424e+00, 3.238113e+00, 5.414722e-02, 4.042938e-01, 1.903040e-01, 6.816790e-01 }   //20-30
+    ,{ 1.745154e-01, 2.234489e+00, 7.962225e-02, 2.007075e+00, 2.925536e+00, 5.499138e-02, 3.958957e-01, 1.780793e+00, 4.852930e-01 }   //30-40
+    ,{ 2.384302e-01, 1.578935e+00, 9.234643e-02, 4.328926e+00, 1.020669e+01, 1.001340e-01, 2.433905e+00, 2.966673e+00, 2.646807e-01 }   //40-50
+    ,{ 0.0000000000, 0.0000000000, 0.0000000000,-1.0000000000, 1.0000000000, 0.0000000000, 0.0000000000, 0.0000000000, 2.0000000000 }   //
+    ,{ 0.0000000000, 0.0000000000, 0.0000000000,-1.0000000000, 1.0000000000, 0.0000000000, 0.0000000000, 0.0000000000, 2.0000000000 }   //pp no v2
+  },{
+    // charged kaon
+    {  0.0000000000, 0.0000000000, 0.0000000000,-1.0000000000, 1.0000000000, 0.0000000000, 0.0000000000, 0.0000000000, 2.0000000000 }  //
+    ,{ 0.0000000000, 0.0000000000, 0.0000000000,-1.0000000000, 1.0000000000, 0.0000000000, 0.0000000000, 0.0000000000, 2.0000000000 }  //
+    ,{ 0.0000000000, 0.0000000000, 0.0000000000,-1.0000000000, 1.0000000000, 0.0000000000, 0.0000000000, 0.0000000000, 2.0000000000 }  //
+    ,{ 1.432314e-01, 1.607297e+00, 2.296868e-01, 1.821156e+00, 2.951652e+00,-2.201385e-01, 1.110419e-01, 8.228701e+00, 4.469488e-01 }  //20-30
+    ,{ 1.691586e-01, 1.617165e+00, 2.159350e-01, 1.649338e+00, 2.607635e+00,-9.536279e+00, 3.860086e-01, 1.802996e+01, 2.509780e-01 }  //30-40
+    ,{ 1.733831e-01, 1.712705e+00, 1.935862e-01, 1.745523e+00, 2.845436e+00,-5.772356e-01, 6.861616e-02, 5.292878e+00, 9.058815e-01 }  //40-50
+    ,{ 0.0000000000, 0.0000000000, 0.0000000000,-1.0000000000, 1.0000000000, 0.0000000000, 0.0000000000, 0.0000000000, 2.0000000000 }  //
+    ,{ 0.0000000000, 0.0000000000, 0.0000000000,-1.0000000000, 1.0000000000, 0.0000000000, 0.0000000000, 0.0000000000, 2.0000000000 }  //pp no v2
+  } };
+
 //==========================================================================
 //
 //              Definition of Particle Distributions
@@ -55,17 +99,23 @@ Double_t AliGenEMlib::PtPizero( const Double_t *px, const Double_t */*dummy*/ )
 {
 // Generate pizero pT distribution from modified Hagedorn parameterization
 // taken from fit to unidentified hadrons in pp at 7 TeV
-  const Double_t kc=0.000565;
-  const Double_t kp0=0.2472;
-  const Double_t kp1=4.354;
-  const Double_t kn=7.007;
+  // const Double_t kc=0.000565;
+  // const Double_t kp0=0.2472;
+  // const Double_t kp1=4.354;
+  // const Double_t kn=7.007;
   Double_t invYield;
   Double_t x=*px;
 
-  invYield = kc/TMath::Power(kp0+x/kp1,kn);
+  // invYield = kc/TMath::Power(kp0+x/kp1,kn);
 
-  return invYield*(2*TMath::Pi()*x);
+  if(!x)return 0;
+  const Double_t *c=fpTparam[fCentbin];
    
+  invYield = CrossOverLc(c[5],c[4],x)*c[0]/TMath::Power(c[3]+x/c[1],c[2]);
+  invYield +=CrossOverRc(c[5],c[4],x)*c[6]/TMath::Power(c[9]+x/c[7],c[8])*CrossOverLc(c[11],c[10],x);
+  invYield +=CrossOverRc(c[11],c[10],x)*c[12]/TMath::Power(x,c[13]);
+
+  return invYield*(2*TMath::Pi()*x);
 }
 
 Double_t AliGenEMlib::YPizero( const Double_t *py, const Double_t */*dummy*/ )
@@ -74,6 +124,11 @@ Double_t AliGenEMlib::YPizero( const Double_t *py, const Double_t */*dummy*/ )
 
 }
 
+Double_t AliGenEMlib::V2Pizero( const Double_t *px, const Double_t */*dummy*/ )
+{
+  return V2Param(px,fv2param[0][fCentbin]);
+}
+
 //--------------------------------------------------------------------------
 //
 //                              Eta
@@ -94,7 +149,11 @@ Double_t AliGenEMlib::PtEta( const Double_t *px, const Double_t */*dummy*/ )
 Double_t AliGenEMlib::YEta( const Double_t *py, const Double_t */*dummy*/ )
 {
   return YFlat(*py);
+}
 
+Double_t AliGenEMlib::V2Eta( const Double_t *px, const Double_t */*dummy*/ )
+{
+  return V2Param(px,fv2param[1][fCentbin]);
 }
 
 //--------------------------------------------------------------------------
@@ -233,6 +292,20 @@ Double_t AliGenEMlib::YFlat(Double_t /*y*/)
   return norm*(pt/scaledPt)*scaledYield;
 }
 
+Double_t AliGenEMlib::V2Param(const Double_t *px, const Double_t *par)
+{
+  // Very general parametrization of the v2
+
+  return TMath::Max(CrossOverLc(par[4],par[3],px[0])*(2*par[0]/(1+TMath::Exp(par[1]*(par[2]-px[0])))-par[0])+CrossOverRc(par[4],par[3],px[0])*((par[8]-par[5])/(1+TMath::Exp(par[6]*(px[0]-par[7])))+par[5]),0.0);
+}
+
+Double_t AliGenEMlib::V2Flat(const Double_t */*px*/, const Double_t */*param*/)
+{
+  // Flat v2
+
+  return 1.0;
+}
+
 //==========================================================================
 //
 //                     Set Getters 
@@ -345,14 +418,36 @@ GenFuncIp AliGenEMlib::GetIp(Int_t param, const char * tname) const
     return func;
 }
 
+GenFunc AliGenEMlib::GetV2(Int_t param, const char * tname) const
+{
+  // Return pointer to v2-parameterisation
+  GenFunc func=0;
+  TString sname(tname);
 
+  switch (param) 
+    {
+    case kPizero:
+      func=V2Pizero;
+      break;
+    case kEta:
+      func=V2Eta;
+      break;
+    case kRho:
+      func=V2Pizero;
+      break;
+    case kOmega:
+      func=V2Pizero;
+      break;
+    case kEtaprime:
+      func=V2Pizero;
+      break;
+    case kPhi:
+      func=V2Pizero;
+      break;
 
-
-
-
-
-
-
-
-
-
+    default:
+      func=0;
+      printf("<AliGenEMlib::GetV2> unknown parametrisation\n");
+    }
+  return func;
+}
index 76ae20b..de51364 100644 (file)
 class TRandom;
 
 class AliGenEMlib :public AliGenLib {
- public:
+public:
     GenFunc   GetPt(Int_t param, const char * tname=0) const;
     GenFunc   GetY(Int_t param, const char * tname=0) const;
     GenFuncIp GetIp(Int_t param, const char * tname=0) const;    
+  GenFunc   GetV2(Int_t param, const char * tname=0) const;
 
-    enum constants{kPizero, kEta, kRho, kOmega, kEtaprime, kPhi};
+  typedef enum {k0005=0, k0510, k1020, k2030, k3040, k4050, k5060, kpp} Centbins_t;
+  enum particles{kPizero, kEta, kRho, kOmega, kEtaprime, kPhi};
 
- private:
+private:
 
-// Pizero
+  // Pizero
     static Int_t    IpPizero(TRandom *ran);
     static Double_t PtPizero( const Double_t *px, const Double_t *dummy );
     static Double_t YPizero(const Double_t *py, const Double_t *dummy);
-    static Double_t MtScal(Double_t pt, Int_t np);
+  static Double_t V2Pizero(const Double_t *px, const Double_t *dummy);
 
-// Eta
+  // Eta
     static Int_t    IpEta(TRandom *ran);
     static Double_t PtEta( const Double_t *px, const Double_t *dummy );
     static Double_t YEta(const Double_t *py, const Double_t *dummy);
+  static Double_t V2Eta(const Double_t *px, const Double_t *dummy);
 
-// Rho
+  // Rho
     static Int_t    IpRho(TRandom *ran);
     static Double_t PtRho( const Double_t *px, const Double_t *dummy );
     static Double_t YRho(const Double_t *py, const Double_t *dummy);
 
-// Omega
+  // Omega
     static Int_t    IpOmega(TRandom *ran);
     static Double_t PtOmega( const Double_t *px, const Double_t *dummy );
     static Double_t YOmega(const Double_t *py, const Double_t *dummy);
 
-// Etaprime
+  // Etaprime
     static Int_t    IpEtaprime(TRandom *ran);
     static Double_t PtEtaprime( const Double_t *px, const Double_t *dummy );
     static Double_t YEtaprime(const Double_t *py, const Double_t *dummy);
 
-// Phi
+  // Phi
     static Int_t    IpPhi(TRandom *ran);
     static Double_t PtPhi( const Double_t *px, const Double_t *dummy );
     static Double_t YPhi(const Double_t *py, const Double_t *dummy);
 
-// General
+  // General
     static Double_t PtFlat(const Double_t *px, const Double_t *dummy);
     static Double_t YFlat(Double_t y);
+  static Double_t MtScal(Double_t pt, Int_t np);
+  static Double_t V2Param(const Double_t *px, const Double_t *param);
+  static Double_t V2Flat(const Double_t *px, const Double_t *param);
 
+  static Double_t CrossOverLc(const double a, const double b, const double x);
+  static Double_t CrossOverRc(const double a, const double b, const double x);
+
+  static const Double_t fv2param[2][8][9];
+  static const Double_t fpTparam[8][14];
+  static Centbins_t fCentbin;
 
   ClassDef(AliGenEMlib,0)
 };
index 306ab15..7b97b67 100644 (file)
@@ -12,21 +12,16 @@ class TRandom;
 class AliGenLib :
   public TObject
 {
- public:
-//
+public:
+  //
     virtual ~AliGenLib(){}
     typedef Double_t (*GenFunc)  (const Double_t *, const Double_t *);
     typedef Int_t    (*GenFuncIp)(TRandom *);    
     virtual GenFunc   GetPt(Int_t param, const char *tname) const   = 0;
     virtual GenFunc   GetY (Int_t param, const char *tname) const   = 0;
     virtual GenFuncIp GetIp(Int_t param, const char *tname) const   = 0;    
+  virtual GenFunc   GetV2(Int_t, const char *) const { return NoV2; }
+  static  Double_t  NoV2(const Double_t *, const Double_t *) { return 0; }
     ClassDef(AliGenLib,0) // Library providing y and pT parameterisations
 };
 #endif
-
-
-
-
-
-
-
index 71847f8..1089601 100644 (file)
@@ -44,19 +44,22 @@ ClassImp(AliGenParam)
 
 //------------------------------------------------------------
 
-  //Begin_Html
-  /*
+//Begin_Html
+/*
     <img src="picts/AliGenParam.gif">
-  */
-  //End_Html
+*/
+//End_Html
 
 //____________________________________________________________
-    AliGenParam::AliGenParam():
-       fPtParaFunc(0),
+AliGenParam::AliGenParam()
+: fPtParaFunc(0),
        fYParaFunc(0),
        fIpParaFunc(0),
+  fV2ParaFunc(0),
        fPtPara(0),
        fYPara(0),
+  fV2Para(0),
+  fdNdPhi(0),
        fParam(0),
        fdNdy0(0.),
        fYWgt(0.),
@@ -67,7 +70,7 @@ ClassImp(AliGenParam)
        fSelectAll(kFALSE),
        fDecayer(0)
 {
-// Default constructor
+  // Default constructor
 }
 //____________________________________________________________
 AliGenParam::AliGenParam(Int_t npart, const AliGenLib * Library,  Int_t param, const char* tname)
@@ -75,8 +78,11 @@ AliGenParam::AliGenParam(Int_t npart, const AliGenLib * Library,  Int_t param, c
      fPtParaFunc(Library->GetPt(param, tname)),
      fYParaFunc (Library->GetY (param, tname)),
      fIpParaFunc(Library->GetIp(param, tname)),
+   fV2ParaFunc(Library->GetV2(param, tname)),
      fPtPara(0),
      fYPara(0),
+   fV2Para(0),
+   fdNdPhi(0),
      fParam(param),
      fdNdy0(0.),
      fYWgt(0.),
@@ -87,7 +93,7 @@ AliGenParam::AliGenParam(Int_t npart, const AliGenLib * Library,  Int_t param, c
      fSelectAll(kFALSE),
      fDecayer(0)
 {
-// Constructor using number of particles parameterisation id and library
+  // Constructor using number of particles parameterisation id and library
     fName = "Param";
     fTitle= "Particle Generator using pT and y parameterisation";
     fAnalog = kAnalog;
@@ -99,8 +105,11 @@ AliGenParam::AliGenParam(Int_t npart, Int_t param, const char* tname, const char
     fPtParaFunc(0),
     fYParaFunc (0),
     fIpParaFunc(0),
+  fV2ParaFunc(0),
     fPtPara(0),
     fYPara(0),
+  fV2Para(0),
+  fdNdPhi(0),
     fParam(param),
     fdNdy0(0.),
     fYWgt(0.),
@@ -111,8 +120,8 @@ AliGenParam::AliGenParam(Int_t npart, Int_t param, const char* tname, const char
     fSelectAll(kFALSE),
     fDecayer(0)
 {
-// Constructor using parameterisation id and number of particles
-//
+  // Constructor using parameterisation id and number of particles
+  //
     fName = name;
     fTitle= "Particle Generator using pT and y parameterisation";
       
@@ -120,6 +129,7 @@ AliGenParam::AliGenParam(Int_t npart, Int_t param, const char* tname, const char
     fPtParaFunc = pLibrary->GetPt(param, tname);
     fYParaFunc  = pLibrary->GetY (param, tname);
     fIpParaFunc = pLibrary->GetIp(param, tname);
+  fV2ParaFunc = pLibrary->GetV2(param, tname);
     
     fAnalog = kAnalog;
     fChildSelect.Set(5);
@@ -136,14 +146,18 @@ AliGenParam::AliGenParam(Int_t npart, Int_t param, const char* tname, const char
 AliGenParam::AliGenParam(Int_t npart, Int_t param,
                          Double_t (*PtPara) (const Double_t*, const Double_t*),
                          Double_t (*YPara ) (const Double_t* ,const Double_t*),
+                         Double_t (*V2Para) (const Double_t* ,const Double_t*),
                         Int_t    (*IpPara) (TRandom *))                 
     :AliGenMC(npart),
      
      fPtParaFunc(PtPara),
      fYParaFunc(YPara),
      fIpParaFunc(IpPara),
+   fV2ParaFunc(V2Para),
      fPtPara(0),
      fYPara(0),
+   fV2Para(0),
+   fdNdPhi(0),
      fParam(param),
      fdNdy0(0.),
      fYWgt(0.),
@@ -154,8 +168,8 @@ AliGenParam::AliGenParam(Int_t npart, Int_t param,
      fSelectAll(kFALSE),
      fDecayer(0)
 {
-// Constructor
-// Gines Martinez 1/10/99 
+  // Constructor
+  // Gines Martinez 1/10/99 
     fName = "Param";
     fTitle= "Particle Generator using pT and y parameterisation";
 
@@ -173,15 +187,17 @@ AliGenParam::AliGenParam(Int_t npart, Int_t param,
 //____________________________________________________________
 AliGenParam::~AliGenParam()
 {
-// Destructor
+  // Destructor
     delete  fPtPara;
     delete  fYPara;
+  delete  fV2Para;
+  delete  fdNdPhi;
 }
 
 //____________________________________________________________
 void AliGenParam::Init()
 {
-// Initialisation
+  // Initialisation
 
     if (gMC) fDecayer = gMC->GetDecayer();
   //Begin_Html
@@ -195,7 +211,7 @@ void AliGenParam::Init()
     if (fPtPara) fPtPara->Delete();
     fPtPara = new TF1(name, fPtParaFunc, fPtMin, fPtMax,0);
     gROOT->GetListOfFunctions()->Remove(fPtPara);
-//  Set representation precision to 10 MeV
+  //  Set representation precision to 10 MeV
     Int_t npx= Int_t((fPtMax - fPtMin) / fDeltaPt);
     
     fPtPara->SetNpx(npx);
@@ -205,24 +221,37 @@ void AliGenParam::Init()
     fYPara  = new TF1(name, fYParaFunc, fYMin, fYMax, 0);
     gROOT->GetListOfFunctions()->Remove(fYPara);
 
+  snprintf(name, 256, "v2-parameterisation for %s", GetName());
+  if (fV2Para) fV2Para->Delete();
+  fV2Para  = new TF1(name, fV2ParaFunc, fPtMin, fPtMax, 0);
+  // fV2Para  = new TF1(name, "2*[0]/(1+TMath::Exp([1]*([2]-x)))-[0]", fPtMin, fPtMax);
+  // fV2Para->SetParameter(0, 0.236910);
+  // fV2Para->SetParameter(1, 1.71122);
+  // fV2Para->SetParameter(2, 0.0827617);
+  //gROOT->GetListOfFunctions()->Remove(fV2Para);  //TR: necessary?
+    
+  snprintf(name, 256, "dNdPhi for %s", GetName());
+  if (fdNdPhi) fdNdPhi->Delete();
+  fdNdPhi  = new TF1(name, "1+2*[0]*TMath::Cos(2*(x-[1]))", fPhiMin, fPhiMax);
+  //gROOT->GetListOfFunctions()->Remove(fdNdPhi);  //TR: necessary?
     
     snprintf(name, 256, "pt-for-%s", GetName());
     TF1 ptPara(name ,fPtParaFunc, 0, 15, 0);
     snprintf(name, 256, "y-for-%s", GetName());
     TF1 yPara(name, fYParaFunc, -6, 6, 0);
 
-//
-// dN/dy| y=0
+  //
+  // dN/dy| y=0
     Double_t y1=0;
     Double_t y2=0;
     
     fdNdy0=fYParaFunc(&y1,&y2);
-//
-// Integral over generation region
+  //
+  // Integral over generation region
     Float_t intYS  = yPara.Integral(fYMin, fYMax,(Double_t*) 0x0,1.e-6);
     Float_t intPt0 = ptPara.Integral(0,15,(Double_t *) 0x0,1.e-6);
     Float_t intPtS = ptPara.Integral(fPtMin,fPtMax,(Double_t*) 0x0,1.e-6);
-    Float_t phiWgt=(fPhiMax-fPhiMin)/2./TMath::Pi();
+  Float_t phiWgt=(fPhiMax-fPhiMin)/2./TMath::Pi();    //TR: should probably be done differently in case of anisotropic phi...
     if (fAnalog == kAnalog) {
        fYWgt  = intYS/fdNdy0;
        fPtWgt = intPtS/intPt0;
@@ -232,34 +261,34 @@ void AliGenParam::Init()
        fPtWgt = (fPtMax-fPtMin)/intPt0;
        fParentWeight = fYWgt*fPtWgt*phiWgt/fNpart;
     }
-//
-// particle decay related initialization
+  //
+  // particle decay related initialization
     fDecayer->SetForceDecay(fForceDecay);
     fDecayer->Init();
 
-//
+  //
     AliGenMC::Init();
 }
 
 //____________________________________________________________
 void AliGenParam::Generate()
 {
-//
-// Generate 'npart' of light and heavy mesons (J/Psi, upsilon or phi, Pion,
-// Kaons, Etas, Omegas) and Baryons (proton, antiprotons, neutrons and 
-// antineutrons in the the desired theta, phi and momentum windows; 
-// Gaussian smearing on the vertex is done if selected. 
-// The decay of heavy mesons is done using lujet, 
-//    and the childern particle are tracked by GEANT
-// However, light mesons are directly tracked by GEANT 
-// setting fForceDecay = nodecay (SetForceDecay(nodecay)) 
-//
-//
-//  Reinitialize decayer
+  //
+  // Generate 'npart' of light and heavy mesons (J/Psi, upsilon or phi, Pion,
+  // Kaons, Etas, Omegas) and Baryons (proton, antiprotons, neutrons and 
+  // antineutrons in the the desired theta, phi and momentum windows; 
+  // Gaussian smearing on the vertex is done if selected. 
+  // The decay of heavy mesons is done using lujet, 
+  //    and the childern particle are tracked by GEANT
+  // However, light mesons are directly tracked by GEANT 
+  // setting fForceDecay = nodecay (SetForceDecay(nodecay)) 
+  //
+  //
+  //  Reinitialize decayer
   fDecayer->SetForceDecay(fForceDecay);
   fDecayer->Init();
 
-//
+  //
   Float_t polar[3]= {0,0,0};  // Polarisation of the parent particle (for GEANT tracking)
   Float_t origin0[3];         // Origin of the generated parent particle (for GEANT tracking)
   Float_t time0;              // Time0 of the generated parent particle
@@ -279,7 +308,7 @@ void AliGenParam::Generate()
   //
   Float_t random[6];
  
-// Calculating vertex position per event
+  // Calculating vertex position per event
   for (j=0;j<3;j++) origin0[j]=fOrigin[j];
   time0 = fTimeOrigin;
   if(fVertexSmear==kPerEvent) {
@@ -290,27 +319,24 @@ void AliGenParam::Generate()
   
   Int_t ipa=0;
   
-// Generating fNpart particles
+  // Generating fNpart particles
   fNprimaries = 0;
   
   while (ipa<fNpart) {
       while(1) {
-//
-// particle type 
+      //
+      // particle type 
          Int_t iPart = fIpParaFunc(fRandom);
          fChildWeight=(fDecayer->GetPartialBranchingRatio(iPart))*fParentWeight;          
          TParticlePDG *particle = pDataBase->GetParticle(iPart);
          Float_t am = particle->Mass();
          
          Rndm(random,2);
-//
-// phi
-         phi=fPhiMin+random[0]*(fPhiMax-fPhiMin);
-//
-// y
+      //
+      // y
          ty = TMath::TanH(fYPara->GetRandom());
-//
-// pT
+      //
+      // pT
          if (fAnalog == kAnalog) {
              pt=fPtPara->GetRandom();
              wgtp=fParentWeight;
@@ -327,15 +353,25 @@ void AliGenParam::Generate()
              Fatal("AliGenParam", 
                    "Division by 0: Please check you rapidity range !");
          }
+      //
+      // phi
+      // if(!ipa)
+      //       phi=fEvPlane; //align first particle of each event with event plane
+      // else{
+       double v2 = fV2Para->Eval(pt);
+       fdNdPhi->SetParameter(0,v2);
+       fdNdPhi->SetParameter(1,fEvPlane);
+       phi=fdNdPhi->GetRandom();
+      // }
          
          pl=xmt*ty/sqrt((1.-ty)*(1.+ty));
          theta=TMath::ATan2(pt,pl);
-// Cut on theta
+      // Cut on theta
          if(theta<fThetaMin || theta>fThetaMax) continue;
          ptot=TMath::Sqrt(pt*pt+pl*pl);
-// Cut on momentum
+      // Cut on momentum
          if(ptot<fPMin || ptot>fPMax) continue;
-//
+      //
          p[0]=pt*TMath::Cos(phi);
          p[1]=pt*TMath::Sin(phi);
          p[2]=pl;
@@ -352,23 +388,23 @@ void AliGenParam::Generate()
                TMath::Sqrt(-2*TMath::Log(random[1]));
          }
          
-// Looking at fForceDecay : 
-// if fForceDecay != none Primary particle decays using 
-// AliPythia and children are tracked by GEANT
-//
-// if fForceDecay == none Primary particle is tracked by GEANT 
-// (In the latest, make sure that GEANT actually does all the decays you want)   
-//
+      // Looking at fForceDecay : 
+      // if fForceDecay != none Primary particle decays using 
+      // AliPythia and children are tracked by GEANT
+      //
+      // if fForceDecay == none Primary particle is tracked by GEANT 
+      // (In the latest, make sure that GEANT actually does all the decays you want)     
+      //
          Bool_t decayed = kFALSE;
          
 
          if (fForceDecay != kNoDecay) {
-// Using lujet to decay particle
+       // Using lujet to decay particle
              Float_t energy=TMath::Sqrt(ptot*ptot+am*am);
              TLorentzVector pmom(p[0], p[1], p[2], energy);
              fDecayer->Decay(iPart,&pmom);
-//
-// select decay particles
+       //
+       // select decay particles
              Int_t np=fDecayer->ImportParticles(particles);
 
              //  Selecting  GeometryAcceptance for particles fPdgCodeParticleforAcceptanceCut;
@@ -395,7 +431,7 @@ void AliGenParam::Generate()
                      iparticle = (TParticle *) particles->At(i);
                      Int_t kf = iparticle->GetPdgCode();
                      Int_t ks = iparticle->GetStatusCode();
-// flagged particle
+           // flagged particle
 
                      if (pFlag[i] == 1) {
                          ipF = iparticle->GetFirstDaughter();
@@ -404,14 +440,14 @@ void AliGenParam::Generate()
                          continue;
                      }
 
-// flag decay products of particles with long life-time (c tau > .3 mum)                     
+           // flag decay products of particles with long life-time (c tau > .3 mum)                  
                      
                      if (ks != 1) { 
-//                       TParticlePDG *particle = pDataBase->GetParticle(kf);
+             //                          TParticlePDG *particle = pDataBase->GetParticle(kf);
                          
                          Double_t lifeTime = fDecayer->GetLifetime(kf);
-//                       Double_t mass     = particle->Mass();
-//                       Double_t width    = particle->Width();
+             //                          Double_t mass     = particle->Mass();
+             //                          Double_t width    = particle->Width();
                          if (lifeTime > (Double_t) fMaxLifeTime) {
                              ipF = iparticle->GetFirstDaughter();
                              ipL = iparticle->GetLastDaughter();       
@@ -421,8 +457,8 @@ void AliGenParam::Generate()
                              pSelected[i]   = 1;
                          }
                      } // ks==1 ?
-//
-// children
+           //
+           // children
                      
                      if ((ChildSelected(TMath::Abs(kf)) || fForceDecay == kAll || fSelectAll) && trackIt[i])
                      {
@@ -449,8 +485,8 @@ void AliGenParam::Generate()
              Int_t iparent;
              if ((fCutOnChild && ncsel >0) || !fCutOnChild){
                  ipa++;
-//
-// Parent
+         //
+         // Parent
                  
                  
                  PushTrack(0, -1, iPart, p, origin0, polar, time0, kPPrimary, nt, wgtp, ((decayed)? 11 : 1));
@@ -458,9 +494,9 @@ void AliGenParam::Generate()
                  KeepTrack(nt); 
                  fNprimaries++;
                  
-//
-// Decay Products
-//               
+         //
+         // Decay Products
+         //              
                  for (i = 1; i < np; i++) {
                      if (pSelected[i]) {
                          TParticle* iparticle = (TParticle *) particles->At(i);
@@ -518,9 +554,9 @@ void AliGenParam::Generate()
 //____________________________________________________________________________________
 Float_t AliGenParam::GetRelativeArea(Float_t ptMin, Float_t ptMax, Float_t yMin, Float_t yMax, Float_t phiMin, Float_t phiMax)
 {
-//
-// Normalisation for selected kinematic region
-//
+  //
+  // Normalisation for selected kinematic region
+  //
   Float_t ratio =  
     fPtPara->Integral(ptMin,ptMax,(Double_t *)0,1.e-6) / fPtPara->Integral( fPtPara->GetXmin(), fPtPara->GetXmax(),(Double_t *)0,1.e-6) *
     fYPara->Integral(yMin,yMax,(Double_t *)0,1.e-6)/fYPara->Integral(fYPara->GetXmin(),fYPara->GetXmax(),(Double_t *)0,1.e-6)   *
index 37a42ad..29c8f7f 100644 (file)
@@ -24,13 +24,14 @@ typedef enum { kAnalog, kNonAnalog} Weighting_t;
 //-------------------------------------------------------------
 class AliGenParam : public AliGenMC
 {
- public:
+public:
     AliGenParam();
     AliGenParam(Int_t npart, const AliGenLib * Library, Int_t param,   const char*  tname = 0);
     AliGenParam(Int_t npart, Int_t param, const char* tname = 0, const char*  name  = 0);
     AliGenParam(Int_t npart, Int_t param,
                Double_t (*PtPara)(const Double_t*, const Double_t*),
                Double_t (*YPara )(const Double_t*, const Double_t*),
+             Double_t (*V2Para)(const Double_t*, const Double_t*),
                Int_t    (*IpPara)(TRandom*)           );
      
     virtual ~AliGenParam();
@@ -49,12 +50,15 @@ class AliGenParam : public AliGenMC
     TF1 *  GetY() {return fYPara;}
     Float_t GetRelativeArea(Float_t ptMin, Float_t ptMax, Float_t yMin, Float_t yMax, Float_t phiMin, Float_t phiMax);
 
- protected:
+protected:
     Double_t (*fPtParaFunc)(const Double_t*, const Double_t*); //! Pointer to Pt parametrisation function
     Double_t (*fYParaFunc )(const Double_t*, const Double_t*); //! Pointer to Y parametrisation function
     Int_t    (*fIpParaFunc )(TRandom*);    //! Pointer to particle type parametrisation function
+  Double_t (*fV2ParaFunc )(const Double_t*, const Double_t*);//! Pointer to V2 parametrisation function
     TF1* fPtPara;              // Transverse momentum parameterisation
     TF1* fYPara;               // Rapidity parameterisation
+  TF1*        fV2Para;       // v2 parametrization
+  TF1*        fdNdPhi;       // Phi distribution depending on v2
     Int_t       fParam;        // Parameterisation type 
     Float_t     fdNdy0;        // central multiplicity per event
     Float_t     fYWgt;         // Y-weight
@@ -65,7 +69,7 @@ class AliGenParam : public AliGenMC
     Bool_t      fSelectAll;    // Flag for transportation of Background while using SetForceDecay()
     AliDecayer  *fDecayer;     // ! Pointer to pythia object for decays
 
- private:
+private:
     AliGenParam(const AliGenParam &Param);
     AliGenParam & operator=(const AliGenParam & rhs);
 
index b856d83..1d082c5 100644 (file)
@@ -81,6 +81,7 @@ AliGenerator::AliGenerator():
   fVertex(3),
   fTimeOrigin(0.),
   fTime(0.),
+  fEvPlane(0),
   fStack(0),
   fContainer(0),
   fCollisionGeometry(0),
@@ -151,6 +152,7 @@ AliGenerator::AliGenerator(Int_t npart):
   fVertex(3),
   fTimeOrigin(0.),
   fTime(0.),
+  fEvPlane(0),
   fStack(0),
   fContainer(0),
   fCollisionGeometry(0),
index fa2aaf0..52e3ae1 100644 (file)
@@ -79,6 +79,7 @@ class AliGenerator : public TNamed, public AliRndm
     virtual TGenerator* GetMC() const {return fMCEvGen;}
     virtual void AddHeader(AliGenEventHeader* /*header*/) {;}
     virtual void SetContainer(AliGenerator* container) {fContainer = container;}
+    virtual void SetEventPlane(Float_t evPlane) {fEvPlane = evPlane; }
 
   // Getters
 
@@ -150,6 +151,7 @@ class AliGenerator : public TNamed, public AliRndm
     
     Float_t     fTimeOrigin; // Time0 origin in a run or event sample
     Float_t     fTime;       // Event time smeared around time0 origin using sigma vertex
+    Float_t     fEvPlane;    // the event plane 
 
     AliStack*   fStack;         //! Local pointer to stack
     AliGenerator* fContainer;   //! Local pointer to container