Added protection agains TRU data.
[u/mrichter/AliRoot.git] / EPOS / AliGenEpos.cxx
CommitLineData
9ef1c2d9 1/*
879c2333 2 * AliGenEpos.cxx
9ef1c2d9 3 *
4 * ALICE event generator based on EPOS model from Klaus Werner
5 *
6 * Created on: Feb 28, 2009
7 * Author: Piotr Ostrowski, postrow@if.pw.edu.pl
8 */
9
10#include "AliGenEpos.h"
11#include "TEpos.h"
12#include "TParticle.h"
879c2333 13#include "TMath.h"
9ef1c2d9 14#include "AliLog.h"
15#include "AliGenEventHeader.h"
16#include "AliGenEposEventHeader.h"
17
18ClassImp(AliGenEpos)
19
20AliGenEpos::AliGenEpos() : AliGenMC(),
21 fBmin(0),
22 fBmax(10000),
23 fPhiMin(0),
879c2333 24 fPhiMax(TMath::TwoPi()),
9ef1c2d9 25 fFilterModelOutput(kFALSE) {
26 SetMC(new TEpos());
27}
28
29AliGenEpos::AliGenEpos(Int_t npart) : AliGenMC(npart),
30 fBmin(0),
31 fBmax(10000),
32 fPhiMin(0),
879c2333 33 fPhiMax(TMath::TwoPi()),
9ef1c2d9 34 fFilterModelOutput(kFALSE) {
35 SetMC(new TEpos());
36}
37
38void AliGenEpos::Init() {
39 AliGenMC::Init();
40 TEpos *epos = GetTEpos();
41 epos->SetLaproj(this->fZProjectile);
42 epos->SetMaproj(this->fAProjectile);
43 epos->SetLatarg(this->fZTarget);
44 epos->SetMatarg(this->fATarget);
45 epos->SetPhimin(this->fPhiMin);
46 epos->SetPhimax(this->fPhiMax);
47 epos->SetBminim(this->fBmin);
48 epos->SetBmaxim(this->fBmax);
49 epos->SetEcms(this->fEnergyCMS);
06c8cc00 50 epos->SetSplitting(kTRUE);
9ef1c2d9 51 GetTEpos()->Initialize();
52}
53
54void AliGenEpos::Generate() {
55 Float_t polar[3] = {0,0,0};
56 Float_t origin0[3] = {0,0,0};
57 Float_t origin[3] = {0,0,0};
879c2333 58 fNprimaries = 0;
9ef1c2d9 59 Int_t nt = 0; //output parameter for PushTrack
60
61 Vertex();
62 for (int j=0; j < 3; j++) origin0[j] = fVertex[j];
63
64 // Generate one event
65
66 GetTEpos()->GenerateEvent();
67 AliWarning("Generated");
68 GetTEpos()->ImportParticles(&fParticles);
69
70 Int_t np = fParticles.GetEntriesFast();
71 AliWarning(Form("Imported %d particles", np));
72
73 Int_t *idsOnStack = NULL;
74 idsOnStack = new Int_t[np];
75 TParticle *iparticle;
76
77 for (int i = 0; i < np; i++) {
78 iparticle = (TParticle *) fParticles.At(i);
79 //Bool_t isNullEntry = iparticle->GetStatusCode() == 0;
80 //Bool_t isCommentOrUnknown = iparticle->GetStatusCode() > 2;
879c2333 81 Bool_t hasDecayed = iparticle->GetStatusCode() >= 2;
9ef1c2d9 82 Bool_t isFinalState = iparticle->GetStatusCode() == 1;
83 Int_t imo = iparticle->GetFirstMother();
84 Bool_t hasMother = (imo >=0);
85
86
87 if (isFinalState) {
88 origin[0] = iparticle->Vx();
89 origin[1] = iparticle->Vy();
90 origin[2] = iparticle->Vz();
91 //doubled track with freeze out coordinates for femtoscopy
92 PushTrack(0,
93 imo>=0?idsOnStack[imo]:-1,
94 iparticle->GetPdgCode(),
95 iparticle->Px(),iparticle->Py(),iparticle->Pz(),iparticle->Energy(),
96 origin[0], origin[1], origin[2],
97 iparticle->T(),
98 polar[0],polar[1],polar[2],
99 hasMother ? kPDecay:kPNoProcess,nt);
100
101 idsOnStack[i] = nt;
102 fNprimaries++;
103 KeepTrack(nt);
104
105 //real track with smeared vertex
106 origin[0] += origin0[0];
107 origin[1] += origin0[1];
108 origin[2] += origin0[2];
109 PushTrack(1,
110 nt, //doubled track as mother
111 iparticle->GetPdgCode(),
112 iparticle->Px(),iparticle->Py(),iparticle->Pz(),iparticle->Energy(),
113 origin[0], origin[1], origin[2],
114 iparticle->T(),
115 polar[0],polar[1],polar[2],
116 kPDecay,nt);
117 fNprimaries++;
118 KeepTrack(nt);
119 } else if(hasDecayed || !fFilterModelOutput) {
120 // don't track it and don't smear vertex
121 origin[0] = iparticle->Vx();
122 origin[1] = iparticle->Vy();
123 origin[2] = iparticle->Vz();
124 PushTrack(0,
125 imo>=0?idsOnStack[imo]:-1,
126 iparticle->GetPdgCode(),
127 iparticle->Px(),iparticle->Py(),iparticle->Pz(),iparticle->Energy(),
128 origin[0], origin[1], origin[2],
129 iparticle->T(),
130 polar[0],polar[1],polar[2],
131 hasMother ? kPDecay:kPNoProcess,nt);
132 idsOnStack[i] = nt;
133 fNprimaries++;
134 KeepTrack(nt);
135 } else if(fFilterModelOutput) {
136 //filtered internal model objects
137 idsOnStack[i] = -1; // to erase mother field in its dauthers
138 }
139 }
140 SetHighWaterMark(fNprimaries);
141 TArrayF eventVertex;
142 eventVertex.Set(3);
143 eventVertex[0] = origin0[0];
144 eventVertex[1] = origin0[1];
145 eventVertex[2] = origin0[2];
146 // Builds the event header, to be called after each event
147 AliGenEposEventHeader* header = new AliGenEposEventHeader("EPOS");
148
149 header->SetNProduced(fNprimaries);
150 header->SetPrimaryVertex(eventVertex);
151
152 header->SetImpactParameter(GetTEpos()->GetBimevt());
153 header->SetReactionPlaneAngle(GetTEpos()->GetPhievt());
154
155 header->SetHardScatters(0);
156 header->SetParticipants(GetTEpos()->GetNpjevt(), GetTEpos()->GetNtgevt());
157 header->SetCollisions(GetTEpos()->GetKolevt(), 0, 0, GetTEpos()->GetNpjevt());
158 header->SetSpectators(GetTEpos()->GetJpnevt(), GetTEpos()->GetJppevt(), GetTEpos()->GetJtnevt(), GetTEpos()->GetJtpevt());
159
160 // Event Vertex
161 header->SetPrimaryVertex(fVertex);
879c2333 162 header->FillInternalFields(GetTEpos());
9ef1c2d9 163 AddHeader(header);
164 fCollisionGeometry = (AliGenEposEventHeader*) header;
165
166 delete[] idsOnStack;
167}
168
169AliGenEpos::~AliGenEpos() {
170}
171