Revised version of the flow afterburner.

It  can boost Hijing events in 3 different ways
(i.e. fixed, random, or hijing r.p.)
Emanuele Simili

diff --git a/EVGEN/AliGenAfterBurnerFlow.cxx b/EVGEN/AliGenAfterBurnerFlow.cxx
index 3bd16e0f3fe1946108d09d8be81011b89cb98777..902360e65070fc442773ad85c52fe366fe378634 100644 (file)
--- a/EVGEN/AliGenAfterBurnerFlow.cxx
+++ b/EVGEN/AliGenAfterBurnerFlow.cxx
@@ -39,6 +39,13 @@
 #include "AliGenAfterBurnerFlow.h"
 #include "AliGenCocktailAfterBurner.h"
 
+// emanuele ---------------------------------------------------------------(
+#include <TList.h>
+#include "AliCollisionGeometry.h"
+#include "AliGenCocktailEntry.h"
+#include "TRandom.h"
+// emanuele ---------------------------------------------------------------)
+
 ClassImp(AliGenAfterBurnerFlow)
 
 ////////////////////////////////////////////////////////////////////////////////////////////////////
@@ -58,14 +65,20 @@ AliGenAfterBurnerFlow::AliGenAfterBurnerFlow(Float_t reactionPlane) {
   //
   // Standard Construction
   // 
-  // reactionPlane - Reaction Plane Angle in Deg [0-360]
-  //
+  // reactionPlane - Reaction Plane Angle given in Deg [0-360]
+  // but stored and applied in radiants (standard for TParticle & AliCollisionGeometry) 
+
+// emanuele ---------------------------------------------------------------(
+
+  if(reactionPlane == 0)     { Info("AliGenAfterBurnerFlow", "Using a random R.P. Angle event by event ( ! not the same used by Hijing ! ) ") ; }
+  else if(reactionPlane < 0) { Info("AliGenAfterBurnerFlow", "Using the Hijing R.P. Angle event by event ") ; }
+  else if(reactionPlane > 0) { Info("AliGenAfterBurnerFlow", "Using a fixed R.P. Angle ( psi = %d deg.) for every event ", reactionPlane) ; }
+  
+  // it was // if(reactionPlane < 0 || reactionPlane > 360)   Error("AliGenAfterBurnerFlow", "Reaction Plane Angle - %d - out of bounds [0-360]", reactionPlane); //
 
-  if (reactionPlane < 0 || reactionPlane > 360)
-    Error("AliGenAfterBurnerFlow", 
-         "Reaction Plane Angle - %d - ot of bounds [0-360]", reactionPlane);
+// emanuele ---------------------------------------------------------------(
 
-  fReactionPlane = reactionPlane;
+  fReactionPlane = 2 * TMath::Pi() * (reactionPlane/360) ;  // r.p. given in degrees (Radomski's way) but stored and applied in radiants (TParticle's way) 
   fCounter = 0;
 }
 
@@ -323,34 +336,95 @@ void AliGenAfterBurnerFlow::Generate() {
 
   // Get Stack of the first Generator
   gen = (AliGenCocktailAfterBurner *)gAlice->Generator();
-  stack = gen->GetStack(0);
-
-  // Loop over particles
-
-  for (Int_t i=0; i<stack->GetNtrack(); i++) {
-
-    particle = stack->Particle(i);
 
-    particle->Momentum(momentum);
-    pdg = particle->GetPdgCode();
-    phi = particle->Phi();
-
-    // get Pt, Y
-
-    pt = momentum.Pt();
-    y = momentum.Rapidity();
-
-    // Calculate Delta Phi for Directed and Elliptic Flow
-    
-    dPhi = -2 * GetCoefficient(pdg, 1, pt, y) * TMath::Sin( phi - fReactionPlane );
-    dPhi -= GetCoefficient(pdg, 2, pt, y) * TMath::Sin( 2 * (phi - fReactionPlane));
-    
-    // Set new phi 
-    
-    phi += dPhi;
-    momentum.SetPhi(phi);
-    particle->SetMomentum(momentum);
+// emanuele ---------------------------------------------------------------(
+
+  AliGenerator* genHijing = 0 ;
+  AliCollisionGeometry* geom = 0 ;
+  AliGenCocktailEntry* entry = 0 ;
+  TList* fEntries = 0 ;
+
+  TRandom* rand = new TRandom(0) ;
+  Float_t fHow = fReactionPlane ;   // this is a temp. solution not to add a new data member in the .h
+
+  for(Int_t ns=0;ns<gen->GetNumberOfEvents();ns++) 
+  {
+   gen->SetActiveEventNumber(ns) ;
+   stack = gen->GetStack(ns);                                  // it was 0.  
+   fEntries = gen->Entries() ;
+
+   TIter next(fEntries) ;
+   while((entry = (AliGenCocktailEntry*)next())) 
+   {
+    if(fHow == 0) // hijing R.P.
+    {
+     Info("Generate (e)","Using R.P. from HIJING ... ");       
+     genHijing = entry->Generator() ;        //  cout <<" * GENERATOR IS "<< genHijing << "  : " << genHijing->GetName() << endl;
+     if(genHijing->ProvidesCollisionGeometry()) 
+     { 
+      geom = gen->GetCollisionGeometry(ns) ; //  cout << " * GEOMETRY YES * " << endl ;
+      fReactionPlane = geom->ReactionPlaneAngle() ; 
+     }
+     else 
+     {
+      Error("Generate (e)", "NO CollisionGeometry !!!  -  using fixed R.P. angle = 0. ") ; 
+      fReactionPlane = 0. ; 
+     }
+    }
+    else if(fHow < 0) // random R.P.
+    {
+     Info("Generate (e)","Using random R.P.s ... ");       
+     fReactionPlane = 2 * TMath::Pi() * rand->Rndm() ;
+    }
+    else  // if constant R.P. -> do nothing (fReactionPlane already setted)
+    {
+     Info("Generate (e)","Using a fixed R.P. psi = %d rad.",fReactionPlane);       
+    }    
+    cout << " * Reaction Plane Angle (event " << ns << ") = " << fReactionPlane << " rad. ( = " << (360*fReactionPlane/(2*TMath::Pi())) << " deg.) * " << endl ;
+   }
+
+// emanuele ---------------------------------------------------------------)
+
+   // Loop over particles
+   
+   for (Int_t i=0; i<stack->GetNtrack(); i++) 
+   {
+     particle = stack->Particle(i);
+
+     particle->Momentum(momentum);
+     pdg = particle->GetPdgCode();
+     phi = particle->Phi();
+
+     // get Pt, Y
+     
+     pt = momentum.Pt() ; 
+     //y = momentum.Rapidity() ;
+
+// emanuele ---------------------------------------------------------------(
+
+    if(TMath::Abs(momentum.Z()) != TMath::Abs(momentum.T())) { y = momentum.Rapidity() ; }
+    else { y = 0. ; }
+    // cout << " * Lorentz Vector (momentum) = " << momentum.X() << " , "  << momentum.Y() << " , " << momentum.Z() << " , " << momentum.T() << " . * " << endl ;
+    // cout << " *                        pt = " << momentum.Pt() << " . * " << endl ;
+    // cout << " *                         Y = " << y << " . * " << endl ;
+
+// emanuele ---------------------------------------------------------------)
+
+     // Calculate Delta Phi for Directed and Elliptic Flow
+     
+     dPhi = -2 * GetCoefficient(pdg, 1, pt, y) * TMath::Sin( phi - fReactionPlane );
+     dPhi -= GetCoefficient(pdg, 2, pt, y) * TMath::Sin( 2 * (phi - fReactionPlane));
+     
+     // Set new phi      
+     
+     phi += dPhi;
+     momentum.SetPhi(phi);
+     particle->SetMomentum(momentum);
+   }
+
+// emanuele ---------------------------------------------------------------(
   }
+// emanuele ---------------------------------------------------------------)
 
   Info("Generate","Flow After Burner: DONE");  
 }

diff --git a/EVGEN/AliGenAfterBurnerFlow.h b/EVGEN/AliGenAfterBurnerFlow.h
index 52aa0e4d952957146f055fa62218a99ca587019f..f37515ac73ba4fe0f3e21a55c214dba038aa88d5 100644 (file)
--- a/EVGEN/AliGenAfterBurnerFlow.h
+++ b/EVGEN/AliGenAfterBurnerFlow.h
@@ -54,7 +54,7 @@ class AliGenAfterBurnerFlow : public AliGenerator {
   Int_t   fCounter;            // counter
 
 
-  ClassDef(AliGenAfterBurnerFlow,2)
+  ClassDef(AliGenAfterBurnerFlow,3)
 
 };
 

diff --git a/EVGEN/AliGenCocktailAfterBurner.cxx b/EVGEN/AliGenCocktailAfterBurner.cxx
index a01c31c18e6dc38542cc113f0dbd509028772510..3f9db9231d08d07ced6c5d92415d7cd7e29c4c22 100644 (file)
--- a/EVGEN/AliGenCocktailAfterBurner.cxx
+++ b/EVGEN/AliGenCocktailAfterBurner.cxx
@@ -39,6 +39,7 @@
 
 #include "AliGenCocktailAfterBurner.h"
 #include "AliGenCocktailEntry.h"
+#include "AliCollisionGeometry.h"
 #include "AliStack.h"
 #include "AliMC.h"
 
@@ -54,7 +55,8 @@ AliGenCocktailAfterBurner::AliGenCocktailAfterBurner()
        cout<<"AliGenCocktailAfterBurner::AliGenCocktailAfterBurner()"<<endl;
     SetName("AliGenCocktailAfterBurner");
     SetTitle("AliGenCocktailAfterBurner");
-    fInternalStacks =0;
+    fInternalStacks = 0;
+    fCollisionGeometries = 0;
     fActiveStack = 0;
     fCurrentGenerator = 0;
     fCurrentEvent =0;
@@ -86,6 +88,7 @@ AliGenCocktailAfterBurner::~AliGenCocktailAfterBurner()
        delete fInternalStacks;
     }
     if (fAfterBurnerEntries) delete fAfterBurnerEntries; //delete entries
+    delete[] fCollisionGeometries;
   }
 /*********************************************************************/ 
 /*********************************************************************/ 
@@ -133,6 +136,10 @@ void AliGenCocktailAfterBurner::Init()
        fInternalStacks->Delete(); //clean after previous generation cycle
      }
 
+// ANDREAS MORSCH ---------------------------------------------------(
+    if (fCollisionGeometries) delete[] fCollisionGeometries;
+// ANDREAS MORSCH ---------------------------------------------------)
+    
     this->AliGenCocktail::Init(); 
     
     if (gDebug>0) cout<<"AliGenCocktailAfterBurner::Init"<<endl;
@@ -179,7 +186,8 @@ void AliGenCocktailAfterBurner::Generate()
       fCurrentEvent=0;
       Int_t numberOfEvents = gAlice->GetEventsPerRun();
       //Create stacks
-      fInternalStacks = new TObjArray(numberOfEvents + fNBgEvents); //Create array of internal stacks
+      fInternalStacks      = new TObjArray(numberOfEvents + fNBgEvents); //Create array of internal stacks
+      fCollisionGeometries = new AliCollisionGeometry*[numberOfEvents + fNBgEvents]; //Create array of collision geometries
       for(i=0;i<numberOfEvents + fNBgEvents;i++) 
        {       
         stack = new AliStack(10000);
@@ -218,6 +226,11 @@ void AliGenCocktailAfterBurner::Generate()
               }
                 fCurrentGenerator->Generate();
                 entry->SetLast(partArray->GetEntriesFast());
+               
+// ANDREAS MORSCH ---------------------------------------------------(
+               if (fCurrentGenerator->ProvidesCollisionGeometry())  fCollisionGeometries[i] = fCurrentGenerator->CollisionGeometry();
+// ANDREAS MORSCH ---------------------------------------------------)
+               
            }
 /***********************************************/
       }
@@ -287,6 +300,25 @@ AliStack* AliGenCocktailAfterBurner::GetStack(Int_t n) const
     }
     return ((AliStack*) fInternalStacks->At(n) );
 }
+
+/*********************************************************************/ 
+/*********************************************************************/ 
+
+// ANDREAS MORSCH ---------------------------------------------------(
+
+AliCollisionGeometry* AliGenCocktailAfterBurner::GetCollisionGeometry(Int_t n) const
+{
+//Returns the pointer to the N'th stack (event)
+  if( ( n<0 ) || ( n>=GetNumberOfEvents() ) )
+    {
+      Fatal("AliGenCocktailAfterBurner::GetCollisionGeometry","Asked for non existing stack (%d)",n);
+      return 0; 
+    }
+    return fCollisionGeometries[n];
+}
+
+// ANDREAS MORSCH ---------------------------------------------------)
+
 /*********************************************************************/ 
 /*********************************************************************/ 
 
@@ -345,8 +377,12 @@ void AliGenCocktailAfterBurner::SetTracks(Int_t stackno)
       mech = AliGenCocktailAfterBurner::IntToMCProcess(p->GetUniqueID());
       weight = p->GetWeight();
 
-      gAlice->GetMCApp()->PushTrack(done, parent, pdg, px, py, pz, e, vx, vy, vz, tof,
-                       polx, poly, polz, mech, ntr, weight);
+      gAlice->GetMCApp()->PushTrack(done, parent, pdg, px, py, pz, e, vx, vy, vz, tof,polx, poly, polz, mech, ntr, weight);
+
+// ANDREAS MORSCH ---------------------------------------------------(
+      SetHighWaterMark(ntr) ; 
+// ANDREAS MORSCH ---------------------------------------------------)
+
     }
 }
 /*********************************************************************/ 

diff --git a/EVGEN/AliGenCocktailAfterBurner.h b/EVGEN/AliGenCocktailAfterBurner.h
index 06e8e452189a677cf3d25a99deefea34408f4ead..c9868516fa569b9c23222bc5bbf19e20fc6d66a7 100644 (file)
--- a/EVGEN/AliGenCocktailAfterBurner.h
+++ b/EVGEN/AliGenCocktailAfterBurner.h
@@ -18,6 +18,10 @@
 class AliGenCocktailEntry;
 class AliStack;
 
+// ANDREAS MORSCH ------------------------------------- (
+class AliCollisionGeometry;
+// ANDREAS MORSCH ------------------------------------- )
+
 class AliGenCocktailAfterBurner : public  AliGenCocktail
 {
 //container class for other generators
@@ -40,7 +44,11 @@ class AliGenCocktailAfterBurner : public  AliGenCocktail
     
     AliStack*     GetStack(Int_t n) const;
     AliStack*     GetActiveStack() const{return fActiveStack;}
-    
+
+// ANDREAS MORSCH ------------------------------------- (
+    AliCollisionGeometry* GetCollisionGeometry(Int_t n) const;
+// ANDREAS MORSCH ------------------------------------- )
+ 
     AliGenerator* GetCurrentGenerator() const;
     virtual void  SetActiveEventNumber(Int_t actev);
     Int_t         GetActiveEventNumber() const {return fActiveEvent;}
@@ -57,6 +65,11 @@ class AliGenCocktailAfterBurner : public  AliGenCocktail
                                 //   if true just return event to gAlice
                                 //   
     TObjArray *fInternalStacks; //! List of internal stacks
+
+// ANDREAS MORSCH ------------------------------------- (
+    AliCollisionGeometry** fCollisionGeometries; //! List of Collision Geometries
+// ANDREAS MORSCH ------------------------------------- )
+    
     Int_t fCurrentEvent;        //  Number of current event/stack