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