EVGEN/AliGenAfterBurnerFlow.cxx

   1
   2 ////////////////////////////////////////////////////////////////////////////////////////////////////
   3 //
   4 // AliGenAfterBurnerFlow is a After Burner event generator applying flow.
   5 // The generator changes Phi coordinate of the particle momentum.
   6 // Flow (directed and elliptical) can be defined on particle type level
   7 //
   8 // For examples, parameters and testing macros refer to:
   9 // http:/home.cern.ch/radomski
  10 //
  11 // Author:
  12 // Sylwester Radomski,
  13 // GSI, April 2002
  14 //
  15 // S.Radomski@gsi.de
  16 //
  17 ////////////////////////////////////////////////////////////////////////////////////////////////////
  18
  19 #include <Riostream.h>
  20 #include "TParticle.h"
  21 #include "TLorentzVector.h"
  22 #include "AliStack.h"
  23 #include "AliGenAfterBurnerFlow.h"
  24 #include "AliGenCocktailAfterBurner.h"
  25
  26 ClassImp(AliGenAfterBurnerFlow)
  27
  28 ////////////////////////////////////////////////////////////////////////////////////////////////////
  29
  30 AliGenAfterBurnerFlow::AliGenAfterBurnerFlow() {
  31   //
  32   // Deafult Construction
  33   //
  34
  35   fReactionPlane = 0;
  36   fCounter = 0;
  37 }
  38
  39 ////////////////////////////////////////////////////////////////////////////////////////////////////
  40
  41 AliGenAfterBurnerFlow::AliGenAfterBurnerFlow(Float_t reactionPlane) {
  42   //
  43   // Standard Construction
  44   //
  45   // reactionPlane - Reaction Plane Angle in Deg [0-360]
  46   //
  47
  48   if (reactionPlane < 0 || reactionPlane > 360)
  49     Error("AliGenAfterBurnerFlow",
  50           "Reaction Plane Angle - %d - ot of bounds [0-360]", reactionPlane);
  51
  52   fReactionPlane = reactionPlane;
  53   fCounter = 0;
  54 }
  55
  56 ////////////////////////////////////////////////////////////////////////////////////////////////////
  57
  58 AliGenAfterBurnerFlow::~AliGenAfterBurnerFlow() {
  59   // Standard Destructor
  60
  61 }
  62
  63 ////////////////////////////////////////////////////////////////////////////////////////////////////
  64
  65 void AliGenAfterBurnerFlow::SetDirectedSimple(Int_t pdg, Float_t v1) {
  66   //
  67   // Set Directed Flow
  68   // The same directed flow is applied to all specified particles
  69   // independently on transverse momentum or rapidity
  70   //
  71   // PDG - particle type to apply directed flow
  72   //       if (PDG == 0) use as default
  73   //
  74
  75   SetFlowParameters(pdg, 1, 0, v1, 0, 0, 0);
  76 }
  77
  78 ////////////////////////////////////////////////////////////////////////////////////////////////////
  79
  80 void AliGenAfterBurnerFlow::SetDirectedParam
  81 (Int_t pdg, Float_t v11, Float_t v12, Float_t v13, Float_t v14) {
  82   //
  83   // Set Directed Flow
  84   // Directed flow is parameterised as follows
  85   //
  86   // V1(Pt,Y) = (V11 + V12*Pt) * sign(Y) * (V13 + V14 * Y^3)
  87   //
  88   // where sign = 1 for Y > 0 and -1 for Y < 0
  89   //
  90   // Defaults values
  91   // v12 = v14 = 0
  92   // v13 = 1
  93   //
  94   // PDG - particle type to apply directed flow
  95   //       if (PDG == 0) use as default
  96   //
  97
  98   SetFlowParameters(pdg, 1, 1, v11, v12, v13, v14);
  99 }
 100
 101 ////////////////////////////////////////////////////////////////////////////////////////////////////
 102
 103 void AliGenAfterBurnerFlow::SetEllipticSimple(Int_t pdg, Float_t v2) {
 104   //
 105   // Set Elliptic Flow
 106   // The same Elliptic flow is applied to all specified particles
 107   // independently on transverse momentum or rapidity
 108   //
 109   // PDG - particle type to apply directed flow
 110   //       if (PDG == 0) use as default
 111   //
 112   // V2 - flow coefficient
 113   //
 114   // NOTE: for starting playing with FLOW
 115   //       start with this function and values 0.05 - 0.1
 116   //
 117
 118   SetFlowParameters(pdg, 2, 0, v2, 0, 0, 0);
 119 }
 120
 121 ////////////////////////////////////////////////////////////////////////////////////////////////////
 122
 123 void AliGenAfterBurnerFlow::SetEllipticParamPion
 124 (Int_t pdg, Float_t v21, Float_t pTmax, Float_t v22) {
 125   //
 126   // Set Elliptic Flow
 127   //
 128   // Elliptic flow is parametrised to reproduce
 129   // V2 of Pions at RHIC energies and is given by:
 130   //
 131   // V2 = v21 * (pT/pTMax ) * exp (-v22 * y^2)    where pT <= pTmax
 132   //      v21 * exp (-v22 * y^2)                   where pT > pTmax
 133   //
 134   // v21   - value at saturation
 135   // pTmax - saturation transverse momentum
 136   // v22   - rapidity decrising
 137   //
 138
 139   SetFlowParameters(pdg, 2, 1, v21, pTmax, v22, 0);
 140 }
 141
 142 ////////////////////////////////////////////////////////////////////////////////////////////////////
 143
 144 void AliGenAfterBurnerFlow::SetEllipticParamOld
 145 (Int_t pdg, Float_t v21, Float_t v22, Float_t v23) {
 146   //
 147   // Set Elliptic Flow
 148   //
 149   // Elliptic flow is parameterised using
 150   // old MevSim parameterisation
 151   //
 152   // V2 = (V21 + V22 pT^2) * exp (-v22 * y^2)
 153   //
 154
 155   SetFlowParameters(pdg, 2, 2, v21, v22, v23, 0);
 156 }
 157
 158 ////////////////////////////////////////////////////////////////////////////////////////////////////
 159
 160 void AliGenAfterBurnerFlow::SetFlowParameters
 161 (Int_t pdg, Int_t order, Int_t type, Float_t v1, Float_t v2,Float_t v3,Float_t v4) {
 162   //
 163   // private function
 164   //
 165
 166   Int_t index = 0;
 167   Bool_t newEntry = kTRUE;
 168
 169   // Defaults
 170
 171   if (pdg == 0) {
 172     index = fgkN - order;
 173     newEntry = kFALSE;
 174   }
 175
 176   // try to find existing entry
 177   for (Int_t i=0; i<fCounter; i++) {
 178     if (pdg == (Int_t)fParams[i][0] &&
 179         order == (Int_t)fParams[i][1]) {
 180
 181       index = i;
 182       newEntry = kFALSE;
 183     }
 184   }
 185
 186   // check fCounter
 187
 188   if (newEntry && (fCounter > fgkN-3)) {
 189     Error("AliAfterBurnerFlow","Overflow");
 190     return;
 191   }
 192
 193   if (newEntry) {
 194     index = fCounter;
 195     fCounter++;
 196   }
 197
 198   // Set new particle type
 199
 200   fParams[index][0] = pdg;
 201   fParams[index][1] = order;
 202   fParams[index][2] = type;
 203   fParams[index][3] = v1;
 204   fParams[index][4] = v2;
 205   fParams[index][5] = v3;
 206   fParams[index][6] = v4;
 207 }
 208
 209 ////////////////////////////////////////////////////////////////////////////////////////////////////
 210
 211 void AliGenAfterBurnerFlow::Init() {
 212   //
 213   // Standard AliGenerator Initializer
 214   //
 215
 216 }
 217
 218 ////////////////////////////////////////////////////////////////////////////////////////////////////
 219
 220 Float_t AliGenAfterBurnerFlow::GetCoefficient
 221 (Int_t pdg, Int_t n, Float_t Pt, Float_t Y) {
 222   //
 223   // private function
 224   // Return Flow Coefficient for a given particle type flow order
 225   // and particle momentum (Pt, Y)
 226   //
 227
 228   Int_t index = fgkN - n;  // default index
 229   Float_t v = 0;
 230
 231   // try to find specific parametrs
 232
 233   for (Int_t i=0; i<fCounter; i++) {
 234
 235     if ((Int_t)fParams[i][0] == pdg &&
 236         (Int_t)fParams[i][1] == n) {
 237
 238       index = i;
 239       break;
 240     }
 241   }
 242
 243   // calculate v
 244
 245   Int_t type = (Int_t)fParams[index][2];
 246
 247   if ((Int_t)fParams[index][1] == 1) { // Directed
 248
 249     if (type == 0 )
 250       v = fParams[index][3];
 251     else
 252       v = (fParams[index][3] + fParams[index][4] * Pt) * TMath::Sign((Float_t)1.,Y) *
 253         (fParams[index][5] + fParams[index][6] * TMath::Abs(Y*Y*Y) );
 254
 255   } else {  // Elliptic
 256
 257     if (type == 0) v = fParams[index][3];
 258
 259     // Pion parameterisation
 260
 261     if (type == 1) {
 262       if (Pt < fParams[index][4])
 263         v = fParams[index][3] * (Pt / fParams[index][4]) ;
 264       else
 265         v = fParams[index][3];
 266
 267       v *= TMath::Exp( - fParams[index][5] * Y * Y);
 268     }
 269
 270     // Old parameterisation
 271
 272     if (type == 2)
 273       v = (fParams[index][3] + fParams[index][4] * Pt * Pt) *
 274         TMath::Exp( - fParams[index][5] * Y * Y);
 275   }
 276
 277   return v;
 278 }
 279
 280 ////////////////////////////////////////////////////////////////////////////////////////////////////
 281
 282 void AliGenAfterBurnerFlow::Generate() {
 283   //
 284   // AliGenerator generate function doing actual job.
 285   // Algorythm:
 286   //
 287   // 1. loop over particles on the stack
 288   // 2. find direct and elliptical flow coefficients for
 289   //    a particle type ore use defaults
 290   // 3. calculate delta phi
 291   // 4. change phi in orginal particle
 292   //
 293   // Algorythm based on:
 294   // A.M. Poskanzer, S.A. Voloshin
 295   // "Methods of analysisng anisotropic flow in relativistic nuclear collisions"
 296   // PRC 58, 1671 (September 1998)
 297   //
 298
 299   AliGenCocktailAfterBurner *gen;
 300   AliStack *stack;
 301   TParticle *particle;
 302   TLorentzVector momentum;
 303
 304   Int_t pdg;
 305   Float_t phi, dPhi;
 306   Float_t pt, y;
 307
 308   // Get Stack of the first Generator
 309   gen = (AliGenCocktailAfterBurner *)gAlice->Generator();
 310   stack = gen->GetStack(0);
 311
 312   // Loop over particles
 313
 314   for (Int_t i=0; i<stack->GetNtrack(); i++) {
 315
 316     particle = stack->Particle(i);
 317
 318     particle->Momentum(momentum);
 319     pdg = particle->GetPdgCode();
 320     phi = particle->Phi();
 321
 322     // get Pt, Y
 323
 324     pt = momentum.Pt();
 325     y = momentum.Rapidity();
 326
 327     // Calculate Delta Phi for Directed and Elliptic Flow
 328
 329     dPhi = -2 * GetCoefficient(pdg, 1, pt, y) * TMath::Sin( phi - fReactionPlane );
 330     dPhi -= GetCoefficient(pdg, 2, pt, y) * TMath::Sin( 2 * (phi - fReactionPlane));
 331
 332     // Set new phi
 333
 334     phi += dPhi;
 335     momentum.SetPhi(phi);
 336     particle->SetMomentum(momentum);
 337   }
 338
 339   Info("Generate","Flow After Burner: DONE");
 340 }
 341
 342 ////////////////////////////////////////////////////////////////////////////////////////////////////
 343