EVGEN/AliGenCocktailAfterBurner.cxx

   1 /**************************************************************************
   2  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
   3  *                                                                        *
   4  * Author: The ALICE Off-line Project.                                    *
   5  * Contributors are mentioned in the code where appropriate.              *
   6  *                                                                        *
   7  * Permission to use, copy, modify and distribute this software and its   *
   8  * documentation strictly for non-commercial purposes is hereby granted   *
   9  * without fee, provided that the above copyright notice appears in all   *
  10  * copies and that both the copyright notice and this permission notice   *
  11  * appear in the supporting documentation. The authors make no claims     *
  12  * about the suitability of this software for any purpose. It is          *
  13  * provided "as is" without express or implied warranty.                  *
  14  **************************************************************************/
  15
  16 //
  17 // Container class for AliGenerator and AfterBurners
  18 // (which are AliGenerators as well) through recursion.
  19 // The container is itself an AliGenerator a
  20 // what is stored are not the pointers to the generators directly
  21 // but to objects of type
  22 // AliGenCocktailAfterBurner entry.
  23 // The class provides also iterator functionality.
  24 // Author: andreas.morsch@cern.ch and piotr.skowronski@cern.ch
  25 //
  26 // 24.09.2001  Piotr Skowronski
  27 //             debug -> gDebug,
  28 //             fNEvents replaced with gAlice->GetEventsPerRun()
  29 //
  30 #include "AliGenCocktailAfterBurner.h"
  31 #include "AliGenCocktailEntry.h"
  32
  33 #include "AliStack.h"
  34 #include <TObjArray.h>
  35 #include <TList.h>
  36 #include <TParticle.h>
  37 #include <Riostream.h>
  38
  39
  40 ClassImp(AliGenCocktailAfterBurner)
  41 /*********************************************************************/
  42 /*********************************************************************/
  43
  44 AliGenCocktailAfterBurner::AliGenCocktailAfterBurner()
  45 {
  46 // Constructor
  47     if (gDebug > 0)
  48        cout<<"AliGenCocktailAfterBurner::AliGenCocktailAfterBurner()"<<endl;
  49     SetName("AliGenCocktailAfterBurner");
  50     SetTitle("AliGenCocktailAfterBurner");
  51     fInternalStacks =0;
  52     fActiveStack = 0;
  53     fCurrentGenerator = 0;
  54     fCurrentEvent =0;
  55     fAfterBurnerEntries = new TList();
  56     fNAfterBurners = 0;
  57     fGenerationDone = kFALSE;
  58
  59     fActiveEvent = -1;
  60     fNBgEvents = 0;
  61 }
  62 /*********************************************************************/
  63 AliGenCocktailAfterBurner::AliGenCocktailAfterBurner(const AliGenCocktailAfterBurner& in)
  64 {
  65  //cpy ctor
  66 }
  67
  68 /*********************************************************************/
  69
  70 AliGenCocktailAfterBurner::~AliGenCocktailAfterBurner()
  71   {
  72 //destructor
  73
  74     if (fInternalStacks) //delete stacks
  75      {
  76        fInternalStacks->SetOwner();
  77        delete fInternalStacks;
  78     }
  79     if (fAfterBurnerEntries) delete fAfterBurnerEntries; //delete entries
  80   }
  81 /*********************************************************************/
  82 /*********************************************************************/
  83
  84 void AliGenCocktailAfterBurner::
  85 AddAfterBurner(AliGenerator *AfterBurner, char* Name, Float_t RateExp)
  86 {
  87 //
  88 //  Forward parameters to the new AfterBurner
  89
  90     if (gDebug>0)cout<<"AliGenCocktailAfterBurner::AddAfterBurner  Named "<<Name<<endl;
  91
  92     if(TestBit(kPtRange))
  93         AfterBurner->SetPtRange(fPtMin,fPtMax);
  94     if(TestBit(kMomentumRange))
  95         AfterBurner->SetMomentumRange(fPMin,fPMax);
  96
  97     AfterBurner->SetYRange(fYMin,fYMax);
  98     AfterBurner->SetPhiRange(fPhiMin*180/TMath::Pi(),fPhiMax*180/TMath::Pi());
  99     AfterBurner->SetThetaRange(fThetaMin*180/TMath::Pi(),fThetaMax*180/TMath::Pi());
 100     AfterBurner->SetOrigin(fOrigin[0], fOrigin[1], fOrigin[2]);
 101     AfterBurner->SetSigma(fOsigma[0], fOsigma[1], fOsigma[2]);
 102     AfterBurner->SetVertexSmear(fVertexSmear);
 103     AfterBurner->SetTrackingFlag(fTrackIt);
 104 //
 105 //  Add AfterBurner to list
 106
 107     AliGenCocktailEntry *entry =
 108         new AliGenCocktailEntry(AfterBurner, Name, RateExp);
 109     fAfterBurnerEntries->Add(entry);
 110     fNAfterBurners++;
 111 //
 112
 113 }
 114 /*********************************************************************/
 115 /*********************************************************************/
 116
 117 void AliGenCocktailAfterBurner::Init()
 118 {
 119 // Initialisation
 120     fGenerationDone = kFALSE;
 121     if (fInternalStacks) //delete stacks
 122      {
 123        fInternalStacks->SetOwner();
 124        fInternalStacks->Delete(); //clean after previous generation cycle
 125      }
 126
 127     this->AliGenCocktail::Init();
 128
 129     if (gDebug>0) cout<<"AliGenCocktailAfterBurner::Init"<<endl;
 130     TIter next(fAfterBurnerEntries);
 131     AliGenCocktailEntry *entry;
 132     //
 133     // Loop over generators and initialize
 134     while((entry = (AliGenCocktailEntry*)next())) {
 135         entry->Generator()->Init();
 136     }
 137 }
 138 /*********************************************************************/
 139 /*********************************************************************/
 140
 141 void AliGenCocktailAfterBurner::Generate()
 142 {
 143 //
 144 // Generate event
 145 //  Firsts runs each generator for all events
 146 //  than after burners ones for each event
 147 //
 148 //  It generates and processes all events during
 149 //  first call only.
 150 //  In next calls it just returns already generated
 151 //  and processed events to the gAlice
 152
 153     if (gDebug>0)
 154       cout<<"#####################################"<<endl
 155           <<"#AliGenCocktailAfterBurner::Generate#"<<endl
 156           <<"#####################################"<<endl;
 157
 158     Int_t i; //iterator
 159     AliStack * stack;
 160
 161     if (fGenerationDone)
 162     {//if generation is done (in first call)
 163      //just copy particles from the stack to the gAlice
 164       SetTracks(++fCurrentEvent);
 165       cout<<"Returning event "<<fCurrentEvent<<endl;
 166       return;
 167     }
 168     else
 169     { //Here we are in the first call of the method
 170       fCurrentEvent=0;
 171       Int_t numberOfEvents = gAlice->GetEventsPerRun();
 172       //Create stacks
 173       fInternalStacks = new TObjArray(numberOfEvents + fNBgEvents); //Create array of internal stacks
 174       for(i=0;i<numberOfEvents + fNBgEvents;i++)
 175        {
 176         stack = new AliStack(10000);
 177         stack->Reset();
 178         fInternalStacks->Add(stack);
 179        }
 180 /*********************************************************************/
 181       TIter next(fEntries);
 182       AliGenCocktailEntry *entry;
 183       AliGenCocktailEntry *e1;
 184       AliGenCocktailEntry *e2;
 185       TObjArray *partArray;
 186   //
 187   // Loop over generators and generate events
 188       Int_t igen=0;
 189       while((entry = (AliGenCocktailEntry*)next()))
 190       {
 191         igen++;
 192         cout<<"Generator "<<igen<<"  : "<<entry->GetName()<<endl;
 193 /***********************************************/
 194 //First generator for all evenets, than second for all events, etc...
 195         for(i=0;i<numberOfEvents + fNBgEvents;i++)
 196           {
 197             cout<<"                  EVENT "<<i<<endl;
 198             stack = GetStack(i);
 199             partArray = stack->Particles();
 200             fCurrentGenerator = entry->Generator();
 201             fCurrentGenerator->SetStack(stack);
 202             if (igen ==1)
 203               {
 204                 entry->SetFirst(0);
 205               }
 206             else
 207               {
 208                 entry->SetFirst((partArray->GetEntriesFast())+1);
 209               }
 210                 fCurrentGenerator->Generate();
 211                 entry->SetLast(partArray->GetEntriesFast());
 212            }
 213 /***********************************************/
 214       }
 215       next.Reset();
 216       while((entry = (AliGenCocktailEntry*)next()))
 217         {
 218           entry->PrintInfo();
 219         }
 220       for ( entry=FirstGenerator();entry;entry=NextGenerator() )
 221         {
 222           entry->PrintInfo();
 223         }
 224       for (FirstGeneratorPair(e1,e2); (e1&&e2); NextGeneratorPair(e1,e2) )
 225         {
 226           printf("\n -----------------------------");
 227           e1->PrintInfo();
 228           e2->PrintInfo();
 229         }
 230
 231
 232       /***********************************************/
 233       /*******After Burners Processing****************/
 234       /***********************************************/
 235       TIter nextAfterBurner(fAfterBurnerEntries);
 236       AliGenCocktailEntry *afterBurnerEntry;
 237       Int_t iab =0; //number of current after burner / counter
 238
 239       cout<<"\n\nRunning After Burners"<<endl;
 240       while((afterBurnerEntry = (AliGenCocktailEntry*)nextAfterBurner()))
 241         {
 242           cout<<"After Burner "<<iab++<<"  :"<<afterBurnerEntry->GetName()<<endl;
 243           fCurrentGenerator = afterBurnerEntry->Generator();
 244           fCurrentGenerator->Generate();
 245         }
 246       cout<<endl<<"Finished. Processed "<<iab<<" After Burners"<<endl;
 247
 248       /***********************************************/
 249       /***********************************************/
 250       /***********************************************/
 251
 252       fGenerationDone=kTRUE;
 253       SetTracks(0); //copy event 0 to gAlice stack
 254
 255 /*********************************************************************/
 256
 257     }//else generated
 258 }
 259 /*********************************************************************/
 260 /*********************************************************************/
 261
 262 AliGenCocktailAfterBurner& AliGenCocktailAfterBurner::operator=(const  AliGenCocktailAfterBurner& rhs)
 263 {
 264 // Assignment operator
 265     return *this;
 266 }
 267 /*********************************************************************/
 268 /*********************************************************************/
 269
 270 AliStack* AliGenCocktailAfterBurner::GetStack(Int_t n) const
 271 {
 272 //Returns the pointer to the N'th stack (event)
 273   if( ( n<0 ) || ( n>=GetNumberOfEvents() ) )
 274     {
 275       Fatal("AliGenCocktailAfterBurner::GetStack","Asked for non existing stack (%d)",n);
 276       return 0;
 277     }
 278     return ((AliStack*) fInternalStacks->At(n) );
 279 }
 280 /*********************************************************************/
 281 /*********************************************************************/
 282
 283 void AliGenCocktailAfterBurner::SetActiveEventNumber(Int_t actev)
 284 {
 285 //Set Active Events Number and Active Stack
 286 //There is only one active event number
 287 //Made fo convinience of work with AfterBurners (HBT processor)
 288
 289     fActiveEvent = actev;
 290     fActiveStack = GetStack(actev);
 291 }
 292 /*********************************************************************/
 293 /*********************************************************************/
 294
 295 void AliGenCocktailAfterBurner::SetTracks(Int_t stackno)
 296 {
 297 //Method which copies tracks from given stack to the
 298 //gAlice's stack
 299     AliStack* instack = GetStack(stackno);
 300     Int_t done;
 301     Int_t parent;
 302     Int_t pdg;
 303     Double_t px, py, pz, e, vx, vy, vz, tof, polx, poly, polz;
 304     TMCProcess mech;
 305     Int_t ntr;
 306     Float_t weight;
 307     TVector3 pol;
 308
 309     TParticle * p;
 310     Int_t n = instack->GetNtrack();
 311     if (gDebug)
 312     {
 313       cout<<"AliGenCocktailAfterBurner::SetTracks("<<stackno<<"). Number of particles is: "<<n<<"\n";
 314     }
 315
 316     for(Int_t i = 0; i < n; i++)
 317     {
 318
 319       p = instack->Particle(i);
 320       done = !p->TestBit(kDoneBit);
 321       parent = p->GetMother(0);
 322       pdg = p->GetPdgCode();
 323       px = p->Px();
 324       py = p->Py();
 325       pz = p->Pz();
 326       e  = p->Energy();
 327       vx = p->Vx();
 328       vy = p->Vy();
 329       vz = p->Vz();
 330       tof = p->T();
 331       p->GetPolarisation(pol);
 332       polx = pol.X();
 333       poly = pol.Y();
 334       polz = pol.Z();
 335       mech = AliGenCocktailAfterBurner::IntToMCProcess(p->GetUniqueID());
 336       weight = p->GetWeight();
 337
 338       gAlice->SetTrack(done, parent, pdg, px, py, pz, e, vx, vy, vz, tof,
 339                        polx, poly, polz, mech, ntr, weight);
 340     }
 341 }
 342 /*********************************************************************/
 343 /*********************************************************************/
 344
 345 TMCProcess AliGenCocktailAfterBurner::IntToMCProcess(Int_t no)
 346 {
 347  //Mothod used to convert uniqueID (integer) to TMCProcess type
 348     const TMCProcess kMCprocesses[kMaxMCProcess] =
 349     {
 350      kPNoProcess, kPMultipleScattering, kPEnergyLoss, kPMagneticFieldL,
 351      kPDecay, kPPair, kPCompton, kPPhotoelectric, kPBrem, kPDeltaRay,
 352      kPAnnihilation, kPHadronic, kPNoProcess, kPEvaporation, kPNuclearFission,
 353      kPNuclearAbsorption, kPPbarAnnihilation, kPNCapture, kPHElastic,
 354      kPHInhelastic, kPMuonNuclear, kPTOFlimit,kPPhotoFission, kPNoProcess,
 355      kPRayleigh, kPNoProcess, kPNoProcess, kPNoProcess, kPNull, kPStop
 356     };
 357
 358     for (Int_t i = 0;i<kMaxMCProcess;i++)
 359     {
 360       if (kMCprocesses[i] == no)
 361         {
 362           return kMCprocesses[i];
 363         }
 364     }
 365     return kPNoProcess;
 366 }