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