[u/mrichter/AliRoot.git] / EVGEN / AliGenCocktailAfterBurner.cxx

/**************************************************************************
 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
 *                                                                        *
 * Author: The ALICE Off-line Project.                                    *
 * Contributors are mentioned in the code where appropriate.              *
 *                                                                        *
 * Permission to use, copy, modify and distribute this software and its   *
 * documentation strictly for non-commercial purposes is hereby granted   *
 * without fee, provided that the above copyright notice appears in all   *
 * copies and that both the copyright notice and this permission notice   *
 * appear in the supporting documentation. The authors make no claims     *
 * about the suitability of this software for any purpose. It is          *
 * provided "as is" without express or implied warranty.                  *
 **************************************************************************/

/* $Id$ */

// 
// Container class for AliGenerator and AfterBurners 
// (which are AliGenerators as well) through recursion.
// The container is itself an AliGenerator a
// what is stored are not the pointers to the generators directly 
// but to objects of type
// AliGenCocktailAfterBurner entry.   
// The class provides also iterator functionality.  
// Author: andreas.morsch@cern.ch and piotr.skowronski@cern.ch
//
// 24.09.2001  Piotr Skowronski
//             debug -> gDebug,
//             fNEvents replaced with AliRunLoader::GetNumberOfEvents()
//


#include <Riostream.h>

#include <TList.h>
#include <TObjArray.h>
#include <TParticle.h>

#include "AliGenCocktailAfterBurner.h"
#include "AliGenCocktailEntry.h"
#include "AliGenCocktailEventHeader.h"
#include "AliCollisionGeometry.h"
#include "AliStack.h"
#include "AliMC.h"
#include "AliRun.h"

using std::cout;
using std::endl;
ClassImp(AliGenCocktailAfterBurner)
/*********************************************************************/ 
/*********************************************************************/ 

    AliGenCocktailAfterBurner::AliGenCocktailAfterBurner():
	fNAfterBurners(0),
	fAfterBurnerEntries(0),
	fGenerationDone(kFALSE),
	fInternalStacks(0),
	fCollisionGeometries(0),
	fHeaders(0),
	fCurrentEvent(0),
	fActiveStack(0),
	fActiveEvent(-1),
	fCurrentGenerator(0),
	fNBgEvents(0)
{
// Constructor
    if (gDebug > 0) 
	cout<<"AliGenCocktailAfterBurner::AliGenCocktailAfterBurner()"<<endl;
    SetName("AliGenCocktailAfterBurner");
    SetTitle("AliGenCocktailAfterBurner");
}

/*********************************************************************/ 

AliGenCocktailAfterBurner::~AliGenCocktailAfterBurner()
  {
//destructor

    if (fInternalStacks) //delete stacks
     { 
       fInternalStacks->SetOwner();
       delete fInternalStacks;
    }
    if (fAfterBurnerEntries) delete fAfterBurnerEntries; //delete entries
    Int_t numberOfEvents = AliRunLoader::Instance()->GetNumberOfEventsPerRun();
    if (fCollisionGeometries) {
      for (Int_t i = 0; i < (numberOfEvents + fNBgEvents); i++)
	if (fCollisionGeometries[i]) delete fCollisionGeometries[i];
      delete[] fCollisionGeometries;
    }
    if (fHeaders) {
      for (Int_t i = 0; i < (numberOfEvents + fNBgEvents); i++)
	if (fHeaders[i]) delete fHeaders[i];
      delete[] fHeaders;
    }
  }
/*********************************************************************/ 
/*********************************************************************/ 

void AliGenCocktailAfterBurner::
AddAfterBurner(AliGenerator *AfterBurner, const char* Name, Float_t RateExp)
{
//
//  Forward parameters to the new AfterBurner
    
    if (gDebug>0)cout<<"AliGenCocktailAfterBurner::AddAfterBurner  Named "<<Name<<endl;

    if(TestBit(kPtRange) && !(AfterBurner->TestBit(kPtRange)) && !(AfterBurner->TestBit(kMomentumRange))) 
	AfterBurner->SetPtRange(fPtMin,fPtMax);
    if(TestBit(kMomentumRange) && !(AfterBurner->TestBit(kPtRange)) && !(AfterBurner->TestBit(kMomentumRange)))
	AfterBurner->SetMomentumRange(fPMin,fPMax);
    
    if (TestBit(kYRange) && !(AfterBurner->TestBit(kYRange)))
	AfterBurner->SetYRange(fYMin,fYMax);
    if (TestBit(kPhiRange) && !(AfterBurner->TestBit(kPhiRange)))
	AfterBurner->SetPhiRange(fPhiMin*180/TMath::Pi(),fPhiMax*180/TMath::Pi());
    if (TestBit(kThetaRange) && !(AfterBurner->TestBit(kThetaRange)) && !(AfterBurner->TestBit(kEtaRange)))
	AfterBurner->SetThetaRange(fThetaMin*180/TMath::Pi(),fThetaMax*180/TMath::Pi());
    if (!(AfterBurner->TestBit(kVertexRange))) {
	AfterBurner->SetOrigin(fOrigin[0], fOrigin[1], fOrigin[2]);
	AfterBurner->SetSigma(fOsigma[0], fOsigma[1], fOsigma[2]);
	AfterBurner->SetVertexSmear(fVertexSmear);
	AfterBurner->SetVertexSource(kContainer);
	AfterBurner->SetTimeOrigin(fTimeOrigin);
    }
    AfterBurner->SetTrackingFlag(fTrackIt);
    //AfterBurner->SetContainer(this);

//
//  Add AfterBurner to list   
    
    AliGenCocktailEntry *entry = 
	new AliGenCocktailEntry(AfterBurner, Name, RateExp);
    if (!fAfterBurnerEntries) fAfterBurnerEntries = new TList();
    
    fAfterBurnerEntries->Add(entry);
    fNAfterBurners++;
//
    
}
/*********************************************************************/ 
/*********************************************************************/ 

void AliGenCocktailAfterBurner::Init()
{
// Initialisation
  Int_t numberOfEvents = AliRunLoader::Instance()->GetNumberOfEventsPerRun();
    fGenerationDone = kFALSE;
    if (fInternalStacks) //delete stacks
     { 
       fInternalStacks->SetOwner();
       fInternalStacks->Delete(); //clean after previous generation cycle
     }

    if (fCollisionGeometries) {
      for (Int_t i = 0; i < (numberOfEvents + fNBgEvents); i++)
	if (fCollisionGeometries[i]) delete fCollisionGeometries[i];
      delete[] fCollisionGeometries;
    }
    if (fHeaders) {
      for (Int_t i = 0; i < (numberOfEvents + fNBgEvents); i++)
	if (fHeaders[i]) delete fHeaders[i];
      delete[] fHeaders;
    }

    this->AliGenCocktail::Init(); 
    
    if (gDebug>0) cout<<"AliGenCocktailAfterBurner::Init"<<endl;
    TIter next(fAfterBurnerEntries);
    AliGenCocktailEntry *entry;
    //
    // Loop over generators and initialize
    while((entry = (AliGenCocktailEntry*)next())) {
	entry->Generator()->Init();
    }  
}
/*********************************************************************/ 
/*********************************************************************/ 

void AliGenCocktailAfterBurner::Generate()
{
//
// Generate event
//  Firsts runs each generator for all events
//  than after burners ones for each event
//
//  It generates and processes all events during
//  first call only.
//  In next calls it just returns already generated 
//  and processed events to the gAlice

    if (gDebug>0)
      cout<<"#####################################"<<endl
          <<"#AliGenCocktailAfterBurner::Generate#"<<endl
          <<"#####################################"<<endl;
    // Initialize header
    //
    Int_t i; //iterator
    AliStack * stack;
    
    if (fGenerationDone)
    {//if generation is done (in first call) 
     //just copy particles from the stack to the gAlice
      SetTracks(++fCurrentEvent);
      fHeader = fHeaders[fCurrentEvent];
      gAlice->SetGenEventHeader(fHeader); 
      cout<<"Returning event " << fCurrentEvent<<endl;
      return;  
    }
    else
    { //Here we are in the first call of the method
	Int_t numberOfEvents = AliRunLoader::Instance()->GetNumberOfEventsPerRun();
	cout << "Number of events per run" <<  numberOfEvents << endl;
	TArrayF eventVertex;
	eventVertex.Set(3 * (numberOfEvents + fNBgEvents));
	TArrayF eventTime;
	eventTime.Set(numberOfEvents + fNBgEvents);
	fCurrentEvent=0;
      //Create stacks
	fInternalStacks      = new TObjArray(numberOfEvents + fNBgEvents); //Create array of internal stacks
	fCollisionGeometries = new AliCollisionGeometry*[numberOfEvents + fNBgEvents]; //Create array of collision geometries
	fHeaders             = new AliGenCocktailEventHeader*[numberOfEvents + fNBgEvents]; //Create array of headers   

	for(i = 0; i < numberOfEvents + fNBgEvents; i++) 
       {	
	   stack = new AliStack(10000);
	   stack->Reset();
	   fInternalStacks->Add(stack);
	   Vertex();
	   for (Int_t j = 0; j < 3; j++) eventVertex[3 * i +  j] = fVertex[j];
	   eventTime[i] = fTime;
	   fHeaders[i] = new AliGenCocktailEventHeader();
	   fCollisionGeometries[i] = 0;
       }
/*********************************************************************/ 
      TIter next(fEntries);
      AliGenCocktailEntry *entry;
      AliGenCocktailEntry *e1;
      AliGenCocktailEntry *e2;
      const TObjArray *partArray;
  //
  // Loop over generators and generate events
      Int_t igen=0;
      while((entry = (AliGenCocktailEntry*)next())) 
      {
        igen++;
        cout<<"Generator "<<igen<<"  : "<<entry->GetName()<<endl;
/***********************************************/
//First generator for all evenets, than second for all events, etc...
        for(i = 0; i < numberOfEvents + fNBgEvents; i++) 
	{  
	    cout<<"                  EVENT "<<i << endl;
            stack = GetStack(i);
            partArray = stack->Particles();
            fCurrentGenerator = entry->Generator();
            fCurrentGenerator->SetStack(stack);
            if (igen ==1) 
	    {
                entry->SetFirst(0);
	    } 
            else 
	    {
                entry->SetFirst((partArray->GetEntriesFast())+1);
	    }
	    // Set the vertex for the generator
	    Int_t ioff = 3 * i;
	    fCurrentGenerator->SetVertex(eventVertex.At(ioff), eventVertex.At(ioff + 1), eventVertex.At(ioff + 2));
	    fCurrentGenerator->SetTime(eventTime.At(i));
	    fHeader = fHeaders[i];
	    // Set the vertex and time for the cocktail
	    TArrayF v(3);
	    for (Int_t j=0; j<3; j++) v[j] = eventVertex.At(ioff + j);
	    fHeader->SetPrimaryVertex(v);
	    fHeader->SetInteractionTime(eventTime.At(i));
	    // Generate event
	    fCurrentGenerator->Generate();
	    //
	    entry->SetLast(partArray->GetEntriesFast());
	    
	    if (fCurrentGenerator->ProvidesCollisionGeometry())  
	      fCollisionGeometries[i] = 
		new AliCollisionGeometry(*(fCurrentGenerator->CollisionGeometry()));
	} // event loop
/***********************************************/
      } // generator loop
      next.Reset();
      while((entry = (AliGenCocktailEntry*)next())) 
        {
          entry->PrintInfo();
        }
      for ( entry=FirstGenerator();entry;entry=NextGenerator() ) 
        {
          entry->PrintInfo();
        }
      for (FirstGeneratorPair(e1,e2); (e1&&e2); NextGeneratorPair(e1,e2) )
        {
          printf("\n -----------------------------");
          e1->PrintInfo();
          e2->PrintInfo();
        }
	
	
      /***********************************************/
      /*******After Burners Processing****************/
      /***********************************************/
      TIter nextAfterBurner(fAfterBurnerEntries);
      AliGenCocktailEntry *afterBurnerEntry;
      Int_t iab =0; //number of current after burner / counter
      
      cout<<"\n\nRunning After Burners"<<endl;
      while((afterBurnerEntry = (AliGenCocktailEntry*)nextAfterBurner()))
        {
          cout<<"After Burner "<<iab++<<"  :"<<afterBurnerEntry->GetName()<<endl;
          fCurrentGenerator = afterBurnerEntry->Generator();
          fCurrentGenerator->Generate();
        }
      cout<<endl<<"Finished. Processed "<<iab<<" After Burners"<<endl;

      /***********************************************/
      /***********************************************/
      /***********************************************/       
        
      fGenerationDone=kTRUE; 
      SetTracks(0); //copy event 0 to gAlice stack
	
/*********************************************************************/
      // Pass the header to gAlice
      fHeader = fHeaders[0];
      gAlice->SetGenEventHeader(fHeader); 
    } //else generated
}
/*********************************************************************/
/*********************************************************************/ 

AliStack* AliGenCocktailAfterBurner::GetStack(Int_t n) const
{
//Returns the pointer to the N'th stack (event)
  if( ( n<0 ) || ( n >= (GetNumberOfEvents()) ) )
    {
      Fatal("AliGenCocktailAfterBurner::GetStack","Asked for non existing stack (%d)",n);
      return 0; 
    }
    return ((AliStack*) fInternalStacks->At(n) );
}

/*********************************************************************/ 
/*********************************************************************/ 

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];
}

/*********************************************************************/ 
/*********************************************************************/ 

void AliGenCocktailAfterBurner::SetActiveEventNumber(Int_t actev)
{
//Set Active Events Number and Active Stack
//There is only one active event number
//Made fo convinience of work with AfterBurners (HBT processor)

    fActiveEvent = actev;
    fActiveStack = GetStack(actev);
}
/*********************************************************************/ 
/*********************************************************************/ 

void AliGenCocktailAfterBurner::SetTracks(Int_t stackno)
{
//Method which copies tracks from given stack to the
//gAlice's stack
    AliStack* instack = GetStack(stackno);
    Int_t done;
    Int_t parent; 
    Int_t pdg;
    Double_t px, py, pz, e, vx, vy, vz, tof, polx, poly, polz;
    TMCProcess mech;
    Int_t ntr;
    Float_t weight;
    TVector3 pol;
    Int_t is;
    
    TParticle * p;
    Int_t n = instack->GetNtrack();
    if (gDebug) 
    {
      cout<<"AliGenCocktailAfterBurner::SetTracks("<<stackno<<"). Number of particles is: "<<n<<"\n";
    }
    
    for(Int_t i = 0; i < n; i++)
    {
	
      p = instack->Particle(i);
      done = !p->TestBit(kDoneBit);
      parent = p->GetMother(0);
      pdg = p->GetPdgCode();
      px = p->Px();
      py = p->Py();
      pz = p->Pz();
      e  = p->Energy();
      vx = p->Vx();
      vy = p->Vy();
      vz = p->Vz();
      tof = p->T();
      p->GetPolarisation(pol);
      polx = pol.X();
      poly = pol.Y();
      polz = pol.Z();
      mech = AliGenCocktailAfterBurner::IntToMCProcess(p->GetUniqueID());
      weight = p->GetWeight();
      is = p->GetStatusCode();

      gAlice->GetMCApp()->PushTrack(done, parent, pdg, px, py, pz, e, vx, vy, vz, tof,polx, poly, polz, mech, ntr, weight, is);

      SetHighWaterMark(ntr) ; 

    }
}
/*********************************************************************/ 
/*********************************************************************/ 

TMCProcess AliGenCocktailAfterBurner::IntToMCProcess(Int_t no)
{
 //Mothod used to convert uniqueID (integer) to TMCProcess type
    const TMCProcess kMCprocesses[kMaxMCProcess] = 
    {
     kPNoProcess, kPMultipleScattering, kPEnergyLoss, kPMagneticFieldL, 
     kPDecay, kPPair, kPCompton, kPPhotoelectric, kPBrem, kPDeltaRay,
     kPAnnihilation, kPHadronic, kPNoProcess, kPEvaporation, kPNuclearFission,
     kPNuclearAbsorption, kPPbarAnnihilation, kPNCapture, kPHElastic, 
     kPHInhelastic, kPMuonNuclear, kPTOFlimit,kPPhotoFission, kPNoProcess, 
     kPRayleigh, kPNoProcess, kPNoProcess, kPNoProcess, kPNull, kPStop
    };
    
    for (Int_t i = 0;i<kMaxMCProcess;i++)
    {
      if (kMCprocesses[i] == no)
        {
          return kMCprocesses[i];
        }
    } 
    return kPNoProcess;
}
Commit	Line	Data
	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 AliRunLoader::GetNumberOfEvents()
	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 "AliGenCocktailEventHeader.h"
	43	#include "AliCollisionGeometry.h"
	44	#include "AliStack.h"
	45	#include "AliMC.h"
	46	#include "AliRun.h"
	47
	48	using std::cout;
	49	using std::endl;
	50	ClassImp(AliGenCocktailAfterBurner)
	51	/*********************************************************************/
	52	/*********************************************************************/
	53
	54	AliGenCocktailAfterBurner::AliGenCocktailAfterBurner():
	55	fNAfterBurners(0),
	56	fAfterBurnerEntries(0),
	57	fGenerationDone(kFALSE),
	58	fInternalStacks(0),
	59	fCollisionGeometries(0),
	60	fHeaders(0),
	61	fCurrentEvent(0),
	62	fActiveStack(0),
	63	fActiveEvent(-1),
	64	fCurrentGenerator(0),
	65	fNBgEvents(0)
	66	{
	67	// Constructor
	68	if (gDebug > 0)
	69	cout<<"AliGenCocktailAfterBurner::AliGenCocktailAfterBurner()"<<endl;
	70	SetName("AliGenCocktailAfterBurner");
	71	SetTitle("AliGenCocktailAfterBurner");
	72	}
	73
	74	/*********************************************************************/
	75
	76	AliGenCocktailAfterBurner::~AliGenCocktailAfterBurner()
	77	{
	78	//destructor
	79
	80	if (fInternalStacks) //delete stacks
	81	{
	82	fInternalStacks->SetOwner();
	83	delete fInternalStacks;
	84	}
	85	if (fAfterBurnerEntries) delete fAfterBurnerEntries; //delete entries
	86	Int_t numberOfEvents = AliRunLoader::Instance()->GetNumberOfEventsPerRun();
	87	if (fCollisionGeometries) {
	88	for (Int_t i = 0; i < (numberOfEvents + fNBgEvents); i++)
	89	if (fCollisionGeometries[i]) delete fCollisionGeometries[i];
	90	delete[] fCollisionGeometries;
	91	}
	92	if (fHeaders) {
	93	for (Int_t i = 0; i < (numberOfEvents + fNBgEvents); i++)
	94	if (fHeaders[i]) delete fHeaders[i];
	95	delete[] fHeaders;
	96	}
	97	}
	98	/*********************************************************************/
	99	/*********************************************************************/
	100
	101	void AliGenCocktailAfterBurner::
	102	AddAfterBurner(AliGenerator AfterBurner, const char Name, Float_t RateExp)
	103	{
	104	//
	105	// Forward parameters to the new AfterBurner
	106
	107	if (gDebug>0)cout<<"AliGenCocktailAfterBurner::AddAfterBurner Named "<<Name<<endl;
	108
	109	if(TestBit(kPtRange) && !(AfterBurner->TestBit(kPtRange)) && !(AfterBurner->TestBit(kMomentumRange)))
	110	AfterBurner->SetPtRange(fPtMin,fPtMax);
	111	if(TestBit(kMomentumRange) && !(AfterBurner->TestBit(kPtRange)) && !(AfterBurner->TestBit(kMomentumRange)))
	112	AfterBurner->SetMomentumRange(fPMin,fPMax);
	113
	114	if (TestBit(kYRange) && !(AfterBurner->TestBit(kYRange)))
	115	AfterBurner->SetYRange(fYMin,fYMax);
	116	if (TestBit(kPhiRange) && !(AfterBurner->TestBit(kPhiRange)))
	117	AfterBurner->SetPhiRange(fPhiMin180/TMath::Pi(),fPhiMax180/TMath::Pi());
	118	if (TestBit(kThetaRange) && !(AfterBurner->TestBit(kThetaRange)) && !(AfterBurner->TestBit(kEtaRange)))
	119	AfterBurner->SetThetaRange(fThetaMin180/TMath::Pi(),fThetaMax180/TMath::Pi());
	120	if (!(AfterBurner->TestBit(kVertexRange))) {
	121	AfterBurner->SetOrigin(fOrigin[0], fOrigin[1], fOrigin[2]);
	122	AfterBurner->SetSigma(fOsigma[0], fOsigma[1], fOsigma[2]);
	123	AfterBurner->SetVertexSmear(fVertexSmear);
	124	AfterBurner->SetVertexSource(kContainer);
	125	AfterBurner->SetTimeOrigin(fTimeOrigin);
	126	}
	127	AfterBurner->SetTrackingFlag(fTrackIt);
	128	//AfterBurner->SetContainer(this);
	129
	130	//
	131	// Add AfterBurner to list
	132
	133	AliGenCocktailEntry *entry =
	134	new AliGenCocktailEntry(AfterBurner, Name, RateExp);
	135	if (!fAfterBurnerEntries) fAfterBurnerEntries = new TList();
	136
	137	fAfterBurnerEntries->Add(entry);
	138	fNAfterBurners++;
	139	//
	140
	141	}
	142	/*********************************************************************/
	143	/*********************************************************************/
	144
	145	void AliGenCocktailAfterBurner::Init()
	146	{
	147	// Initialisation
	148	Int_t numberOfEvents = AliRunLoader::Instance()->GetNumberOfEventsPerRun();
	149	fGenerationDone = kFALSE;
	150	if (fInternalStacks) //delete stacks
	151	{
	152	fInternalStacks->SetOwner();
	153	fInternalStacks->Delete(); //clean after previous generation cycle
	154	}
	155
	156	if (fCollisionGeometries) {
	157	for (Int_t i = 0; i < (numberOfEvents + fNBgEvents); i++)
	158	if (fCollisionGeometries[i]) delete fCollisionGeometries[i];
	159	delete[] fCollisionGeometries;
	160	}
	161	if (fHeaders) {
	162	for (Int_t i = 0; i < (numberOfEvents + fNBgEvents); i++)
	163	if (fHeaders[i]) delete fHeaders[i];
	164	delete[] fHeaders;
	165	}
	166
	167	this->AliGenCocktail::Init();
	168
	169	if (gDebug>0) cout<<"AliGenCocktailAfterBurner::Init"<<endl;
	170	TIter next(fAfterBurnerEntries);
	171	AliGenCocktailEntry *entry;
	172	//
	173	// Loop over generators and initialize
	174	while((entry = (AliGenCocktailEntry*)next())) {
	175	entry->Generator()->Init();
	176	}
	177	}
	178	/*********************************************************************/
	179	/*********************************************************************/
	180
	181	void AliGenCocktailAfterBurner::Generate()
	182	{
	183	//
	184	// Generate event
	185	// Firsts runs each generator for all events
	186	// than after burners ones for each event
	187	//
	188	// It generates and processes all events during
	189	// first call only.
	190	// In next calls it just returns already generated
	191	// and processed events to the gAlice
	192
	193	if (gDebug>0)
	194	cout<<"#####################################"<<endl
	195	<<"#AliGenCocktailAfterBurner::Generate#"<<endl
	196	<<"#####################################"<<endl;
	197	// Initialize header
	198	//
	199	Int_t i; //iterator
	200	AliStack * stack;
	201
	202	if (fGenerationDone)
	203	{//if generation is done (in first call)
	204	//just copy particles from the stack to the gAlice
	205	SetTracks(++fCurrentEvent);
	206	fHeader = fHeaders[fCurrentEvent];
	207	gAlice->SetGenEventHeader(fHeader);
	208	cout<<"Returning event " << fCurrentEvent<<endl;
	209	return;
	210	}
	211	else
	212	{ //Here we are in the first call of the method
	213	Int_t numberOfEvents = AliRunLoader::Instance()->GetNumberOfEventsPerRun();
	214	cout << "Number of events per run" << numberOfEvents << endl;
	215	TArrayF eventVertex;
	216	eventVertex.Set(3 * (numberOfEvents + fNBgEvents));
	217	TArrayF eventTime;
	218	eventTime.Set(numberOfEvents + fNBgEvents);
	219	fCurrentEvent=0;
	220	//Create stacks
	221	fInternalStacks = new TObjArray(numberOfEvents + fNBgEvents); //Create array of internal stacks
	222	fCollisionGeometries = new AliCollisionGeometry*[numberOfEvents + fNBgEvents]; //Create array of collision geometries
	223	fHeaders = new AliGenCocktailEventHeader*[numberOfEvents + fNBgEvents]; //Create array of headers
	224
	225	for(i = 0; i < numberOfEvents + fNBgEvents; i++)
	226	{
	227	stack = new AliStack(10000);
	228	stack->Reset();
	229	fInternalStacks->Add(stack);
	230	Vertex();
	231	for (Int_t j = 0; j < 3; j++) eventVertex[3 * i + j] = fVertex[j];
	232	eventTime[i] = fTime;
	233	fHeaders[i] = new AliGenCocktailEventHeader();
	234	fCollisionGeometries[i] = 0;
	235	}
	236	/*********************************************************************/
	237	TIter next(fEntries);
	238	AliGenCocktailEntry *entry;
	239	AliGenCocktailEntry *e1;
	240	AliGenCocktailEntry *e2;
	241	const TObjArray *partArray;
	242	//
	243	// Loop over generators and generate events
	244	Int_t igen=0;
	245	while((entry = (AliGenCocktailEntry*)next()))
	246	{
	247	igen++;
	248	cout<<"Generator "<<igen<<" : "<<entry->GetName()<<endl;
	249	/***********************************************/
	250	//First generator for all evenets, than second for all events, etc...
	251	for(i = 0; i < numberOfEvents + fNBgEvents; i++)
	252	{
	253	cout<<" EVENT "<<i << endl;
	254	stack = GetStack(i);
	255	partArray = stack->Particles();
	256	fCurrentGenerator = entry->Generator();
	257	fCurrentGenerator->SetStack(stack);
	258	if (igen ==1)
	259	{
	260	entry->SetFirst(0);
	261	}
	262	else
	263	{
	264	entry->SetFirst((partArray->GetEntriesFast())+1);
	265	}
	266	// Set the vertex for the generator
	267	Int_t ioff = 3 * i;
	268	fCurrentGenerator->SetVertex(eventVertex.At(ioff), eventVertex.At(ioff + 1), eventVertex.At(ioff + 2));
	269	fCurrentGenerator->SetTime(eventTime.At(i));
	270	fHeader = fHeaders[i];
	271	// Set the vertex and time for the cocktail
	272	TArrayF v(3);
	273	for (Int_t j=0; j<3; j++) v[j] = eventVertex.At(ioff + j);
	274	fHeader->SetPrimaryVertex(v);
	275	fHeader->SetInteractionTime(eventTime.At(i));
	276	// Generate event
	277	fCurrentGenerator->Generate();
	278	//
	279	entry->SetLast(partArray->GetEntriesFast());
	280
	281	if (fCurrentGenerator->ProvidesCollisionGeometry())
	282	fCollisionGeometries[i] =
	283	new AliCollisionGeometry(*(fCurrentGenerator->CollisionGeometry()));
	284	} // event loop
	285	/***********************************************/
	286	} // generator loop
	287	next.Reset();
	288	while((entry = (AliGenCocktailEntry*)next()))
	289	{
	290	entry->PrintInfo();
	291	}
	292	for ( entry=FirstGenerator();entry;entry=NextGenerator() )
	293	{
	294	entry->PrintInfo();
	295	}
	296	for (FirstGeneratorPair(e1,e2); (e1&&e2); NextGeneratorPair(e1,e2) )
	297	{
	298	printf("\n -----------------------------");
	299	e1->PrintInfo();
	300	e2->PrintInfo();
	301	}
	302
	303
	304	/***********************************************/
	305	/*****After Burners Processing**************/
	306	/***********************************************/
	307	TIter nextAfterBurner(fAfterBurnerEntries);
	308	AliGenCocktailEntry *afterBurnerEntry;
	309	Int_t iab =0; //number of current after burner / counter
	310
	311	cout<<"\n\nRunning After Burners"<<endl;
	312	while((afterBurnerEntry = (AliGenCocktailEntry*)nextAfterBurner()))
	313	{
	314	cout<<"After Burner "<<iab++<<" :"<<afterBurnerEntry->GetName()<<endl;
	315	fCurrentGenerator = afterBurnerEntry->Generator();
	316	fCurrentGenerator->Generate();
	317	}
	318	cout<<endl<<"Finished. Processed "<<iab<<" After Burners"<<endl;
	319
	320	/***********************************************/
	321	/***********************************************/
	322	/***********************************************/
	323
	324	fGenerationDone=kTRUE;
	325	SetTracks(0); //copy event 0 to gAlice stack
	326
	327	/*********************************************************************/
	328	// Pass the header to gAlice
	329	fHeader = fHeaders[0];
	330	gAlice->SetGenEventHeader(fHeader);
	331	} //else generated
	332	}
	333	/*********************************************************************/
	334	/*********************************************************************/
	335
	336	AliStack* AliGenCocktailAfterBurner::GetStack(Int_t n) const
	337	{
	338	//Returns the pointer to the N'th stack (event)
	339	if( ( n<0 ) \|\| ( n >= (GetNumberOfEvents()) ) )
	340	{
	341	Fatal("AliGenCocktailAfterBurner::GetStack","Asked for non existing stack (%d)",n);
	342	return 0;
	343	}
	344	return ((AliStack*) fInternalStacks->At(n) );
	345	}
	346
	347	/*********************************************************************/
	348	/*********************************************************************/
	349
	350	AliCollisionGeometry* AliGenCocktailAfterBurner::GetCollisionGeometry(Int_t n) const
	351	{
	352	//Returns the pointer to the N'th stack (event)
	353	if( ( n<0 ) \|\| ( n>=GetNumberOfEvents() ) )
	354	{
	355	Fatal("AliGenCocktailAfterBurner::GetCollisionGeometry","Asked for non existing stack (%d)",n);
	356	return 0;
	357	}
	358	return fCollisionGeometries[n];
	359	}
	360
	361	/*********************************************************************/
	362	/*********************************************************************/
	363
	364	void AliGenCocktailAfterBurner::SetActiveEventNumber(Int_t actev)
	365	{
	366	//Set Active Events Number and Active Stack
	367	//There is only one active event number
	368	//Made fo convinience of work with AfterBurners (HBT processor)
	369
	370	fActiveEvent = actev;
	371	fActiveStack = GetStack(actev);
	372	}
	373	/*********************************************************************/
	374	/*********************************************************************/
	375
	376	void AliGenCocktailAfterBurner::SetTracks(Int_t stackno)
	377	{
	378	//Method which copies tracks from given stack to the
	379	//gAlice's stack
	380	AliStack* instack = GetStack(stackno);
	381	Int_t done;
	382	Int_t parent;
	383	Int_t pdg;
	384	Double_t px, py, pz, e, vx, vy, vz, tof, polx, poly, polz;
	385	TMCProcess mech;
	386	Int_t ntr;
	387	Float_t weight;
	388	TVector3 pol;
	389	Int_t is;
	390
	391	TParticle * p;
	392	Int_t n = instack->GetNtrack();
	393	if (gDebug)
	394	{
	395	cout<<"AliGenCocktailAfterBurner::SetTracks("<<stackno<<"). Number of particles is: "<<n<<"\n";
	396	}
	397
	398	for(Int_t i = 0; i < n; i++)
	399	{
	400
	401	p = instack->Particle(i);
	402	done = !p->TestBit(kDoneBit);
	403	parent = p->GetMother(0);
	404	pdg = p->GetPdgCode();
	405	px = p->Px();
	406	py = p->Py();
	407	pz = p->Pz();
	408	e = p->Energy();
	409	vx = p->Vx();
	410	vy = p->Vy();
	411	vz = p->Vz();
	412	tof = p->T();
	413	p->GetPolarisation(pol);
	414	polx = pol.X();
	415	poly = pol.Y();
	416	polz = pol.Z();
	417	mech = AliGenCocktailAfterBurner::IntToMCProcess(p->GetUniqueID());
	418	weight = p->GetWeight();
	419	is = p->GetStatusCode();
	420
	421	gAlice->GetMCApp()->PushTrack(done, parent, pdg, px, py, pz, e, vx, vy, vz, tof,polx, poly, polz, mech, ntr, weight, is);
	422
	423	SetHighWaterMark(ntr) ;
	424
	425	}
	426	}
	427	/*********************************************************************/
	428	/*********************************************************************/
	429
	430	TMCProcess AliGenCocktailAfterBurner::IntToMCProcess(Int_t no)
	431	{
	432	//Mothod used to convert uniqueID (integer) to TMCProcess type
	433	const TMCProcess kMCprocesses[kMaxMCProcess] =
	434	{
	435	kPNoProcess, kPMultipleScattering, kPEnergyLoss, kPMagneticFieldL,
	436	kPDecay, kPPair, kPCompton, kPPhotoelectric, kPBrem, kPDeltaRay,
	437	kPAnnihilation, kPHadronic, kPNoProcess, kPEvaporation, kPNuclearFission,
	438	kPNuclearAbsorption, kPPbarAnnihilation, kPNCapture, kPHElastic,
	439	kPHInhelastic, kPMuonNuclear, kPTOFlimit,kPPhotoFission, kPNoProcess,
	440	kPRayleigh, kPNoProcess, kPNoProcess, kPNoProcess, kPNull, kPStop
	441	};
	442
	443	for (Int_t i = 0;i<kMaxMCProcess;i++)
	444	{
	445	if (kMCprocesses[i] == no)
	446	{
	447	return kMCprocesses[i];
	448	}
	449	}
	450	return kPNoProcess;
	451	}
	452