]>
Commit | Line | Data |
---|---|---|
e40b9538 | 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: AliGenEMCocktail.cxx 40702 2010-04-26 13:09:52Z morsch $ */ | |
17 | ||
4ae1c9f0 | 18 | // Class to create the cocktail for physics with electrons, di-electrons, |
e40b9538 | 19 | // and photons from the decay of the following sources: |
20 | // pizero, eta, rho, omega, etaprime, phi | |
21 | // Kinematic distributions of the sources are taken from AliGenEMlib. | |
22 | // Decay channels can be selected via the method SetDecayMode. | |
23 | // Particles can be generated flat in pT with weights according to the | |
24 | // chosen pT distributions from AliGenEMlib (weighting mode: kNonAnalog), | |
25 | // or they are generated according to the pT distributions themselves | |
26 | // (weighting mode: kAnalog) | |
27 | ||
4ae1c9f0 | 28 | |
e40b9538 | 29 | #include <TObjArray.h> |
30 | #include <TParticle.h> | |
31 | #include <TF1.h> | |
32 | #include <TVirtualMC.h> | |
33 | #include <TPDGCode.h> | |
34 | #include "AliGenCocktailEventHeader.h" | |
35 | ||
36 | #include "AliGenCocktailEntry.h" | |
37 | #include "AliGenEMCocktail.h" | |
38 | #include "AliGenEMlib.h" | |
39 | #include "AliGenBox.h" | |
40 | #include "AliGenParam.h" | |
41 | #include "AliMC.h" | |
42 | #include "AliRun.h" | |
43 | #include "AliStack.h" | |
44 | #include "AliDecayer.h" | |
45 | #include "AliDecayerPythia.h" | |
46 | #include "AliLog.h" | |
47 | #include "AliGenCorrHF.h" | |
48 | ||
49 | ClassImp(AliGenEMCocktail) | |
50 | ||
51 | //________________________________________________________________________ | |
52 | AliGenEMCocktail::AliGenEMCocktail() | |
4ae1c9f0 | 53 | :AliGenCocktail(), |
e40b9538 | 54 | fDecayer(0), |
55 | fDecayMode(kAll), | |
56 | fWeightingMode(kNonAnalog), | |
57 | fNPart(1000), | |
4ae1c9f0 | 58 | fYieldArray(), |
59 | fPtSelect(0), | |
60 | fCentrality(0), | |
61 | fV2Systematic(0), | |
62 | fForceConv(kFALSE), | |
63 | fHeaviestParticle(kGENs) | |
e40b9538 | 64 | { |
4ae1c9f0 | 65 | // Constructor |
e40b9538 | 66 | |
67 | } | |
68 | ||
69 | //_________________________________________________________________________ | |
70 | AliGenEMCocktail::~AliGenEMCocktail() | |
71 | { | |
4ae1c9f0 | 72 | // Destructor |
e40b9538 | 73 | |
74 | } | |
75 | ||
76 | //_________________________________________________________________________ | |
77 | void AliGenEMCocktail::CreateCocktail() | |
78 | { | |
4ae1c9f0 | 79 | // create and add sources to the cocktail |
e40b9538 | 80 | |
81 | fDecayer->SetForceDecay(fDecayMode); | |
82 | fDecayer->ForceDecay(); | |
83 | ||
4ae1c9f0 | 84 | // Set kinematic limits |
e40b9538 | 85 | Double_t ptMin = fPtMin; |
86 | Double_t ptMax = fPtMax; | |
87 | Double_t yMin = fYMin;; | |
88 | Double_t yMax = fYMax;; | |
89 | Double_t phiMin = fPhiMin*180./TMath::Pi(); | |
90 | Double_t phiMax = fPhiMax*180./TMath::Pi(); | |
91 | AliInfo(Form("Ranges pT:%4.1f : %4.1f GeV/c, y:%4.2f : %4.2f, Phi:%5.1f : %5.1f degres",ptMin,ptMax,yMin,yMax,phiMin,phiMax)); | |
92 | AliInfo(Form("the parametrised sources uses the decay mode %d",fDecayMode)); | |
93 | ||
4ae1c9f0 | 94 | //Initialize user selection for Pt Parameterization and centrality: |
95 | AliGenEMlib::SelectParams(fPtSelect,fCentrality,fV2Systematic); | |
96 | ||
97 | // Create and add electron sources to the generator | |
e40b9538 | 98 | |
4ae1c9f0 | 99 | // pizero |
100 | AliGenParam *genpizero=0; | |
e40b9538 | 101 | Char_t namePizero[10]; |
4ae1c9f0 | 102 | snprintf(namePizero,10,"Pizero"); |
e40b9538 | 103 | genpizero = new AliGenParam(fNPart, new AliGenEMlib(), AliGenEMlib::kPizero, "DUMMY"); |
104 | AddSource2Generator(namePizero,genpizero); | |
105 | TF1 *fPtPizero = genpizero->GetPt(); | |
faf00b0d | 106 | #if ROOT_VERSION_CODE >= ROOT_VERSION(5,99,0) |
107 | fYieldArray[kGenPizero] = fPtPizero->Integral(fPtMin,fPtMax,1.e-6); | |
108 | #else | |
4ae1c9f0 | 109 | fYieldArray[kGenPizero] = fPtPizero->Integral(fPtMin,fPtMax,(Double_t *)0,1.e-6); |
faf00b0d | 110 | #endif |
4ae1c9f0 | 111 | |
112 | // eta | |
113 | if(fHeaviestParticle<kGenEta)return; | |
114 | AliGenParam *geneta=0; | |
e40b9538 | 115 | Char_t nameEta[10]; |
4ae1c9f0 | 116 | snprintf(nameEta,10,"Eta"); |
e40b9538 | 117 | geneta = new AliGenParam(fNPart, new AliGenEMlib(), AliGenEMlib::kEta, "DUMMY"); |
118 | AddSource2Generator(nameEta,geneta); | |
119 | TF1 *fPtEta = geneta->GetPt(); | |
4ae1c9f0 | 120 | #if ROOT_VERSION_CODE >= ROOT_VERSION(5,99,0) |
faf00b0d | 121 | fYieldArray[kGenEta] = fPtEta->Integral(fPtMin,fPtMax,1.e-6); |
4ae1c9f0 | 122 | #else |
123 | fYieldArray[kGenEta] = fPtEta->Integral(fPtMin,fPtMax,(Double_t *)0,1.e-6); | |
faf00b0d | 124 | #endif |
4ae1c9f0 | 125 | |
126 | // rho | |
127 | if(fHeaviestParticle<kGenRho)return; | |
128 | AliGenParam *genrho=0; | |
e40b9538 | 129 | Char_t nameRho[10]; |
4ae1c9f0 | 130 | snprintf(nameRho,10,"Rho"); |
e40b9538 | 131 | genrho = new AliGenParam(fNPart, new AliGenEMlib(), AliGenEMlib::kRho, "DUMMY"); |
132 | AddSource2Generator(nameRho,genrho); | |
133 | TF1 *fPtRho = genrho->GetPt(); | |
faf00b0d | 134 | #if ROOT_VERSION_CODE >= ROOT_VERSION(5,99,0) |
135 | fYieldArray[kGenRho] = fPtRho->Integral(fPtMin,fPtMax,1.e-6); | |
136 | #else | |
4ae1c9f0 | 137 | fYieldArray[kGenRho] = fPtRho->Integral(fPtMin,fPtMax,(Double_t *)0,1.e-6); |
faf00b0d | 138 | #endif |
4ae1c9f0 | 139 | |
140 | // omega | |
141 | if(fHeaviestParticle<kGenOmega)return; | |
142 | AliGenParam *genomega=0; | |
e40b9538 | 143 | Char_t nameOmega[10]; |
4ae1c9f0 | 144 | snprintf(nameOmega,10,"Omega"); |
e40b9538 | 145 | genomega = new AliGenParam(fNPart, new AliGenEMlib(), AliGenEMlib::kOmega, "DUMMY"); |
146 | AddSource2Generator(nameOmega,genomega); | |
147 | TF1 *fPtOmega = genomega->GetPt(); | |
faf00b0d | 148 | #if ROOT_VERSION_CODE >= ROOT_VERSION(5,99,0) |
149 | fYieldArray[kGenOmega] = fPtOmega->Integral(fPtMin,fPtMax,1.e-6); | |
150 | #else | |
4ae1c9f0 | 151 | fYieldArray[kGenOmega] = fPtOmega->Integral(fPtMin,fPtMax,(Double_t *)0,1.e-6); |
faf00b0d | 152 | #endif |
4ae1c9f0 | 153 | |
154 | // etaprime | |
155 | if(fHeaviestParticle<kGenEtaprime)return; | |
156 | AliGenParam *genetaprime=0; | |
e40b9538 | 157 | Char_t nameEtaprime[10]; |
4ae1c9f0 | 158 | snprintf(nameEtaprime,10,"Etaprime"); |
e40b9538 | 159 | genetaprime = new AliGenParam(fNPart, new AliGenEMlib(), AliGenEMlib::kEtaprime, "DUMMY"); |
160 | AddSource2Generator(nameEtaprime,genetaprime); | |
161 | TF1 *fPtEtaprime = genetaprime->GetPt(); | |
faf00b0d | 162 | #if ROOT_VERSION_CODE >= ROOT_VERSION(5,99,0) |
163 | fYieldArray[kGenEtaprime] = fPtEtaprime->Integral(fPtMin,fPtMax,1.e-6); | |
164 | #else | |
4ae1c9f0 | 165 | fYieldArray[kGenEtaprime] = fPtEtaprime->Integral(fPtMin,fPtMax,(Double_t *)0,1.e-6); |
faf00b0d | 166 | #endif |
4ae1c9f0 | 167 | |
168 | // phi | |
169 | if(fHeaviestParticle<kGenPhi)return; | |
170 | AliGenParam *genphi=0; | |
e40b9538 | 171 | Char_t namePhi[10]; |
4ae1c9f0 | 172 | snprintf(namePhi,10,"Phi"); |
e40b9538 | 173 | genphi = new AliGenParam(fNPart, new AliGenEMlib(), AliGenEMlib::kPhi, "DUMMY"); |
174 | AddSource2Generator(namePhi,genphi); | |
175 | TF1 *fPtPhi = genphi->GetPt(); | |
faf00b0d | 176 | #if ROOT_VERSION_CODE >= ROOT_VERSION(5,99,0) |
177 | fYieldArray[kGenPhi] = fPtPhi->Integral(fPtMin,fPtMax,1.e-6); | |
178 | #else | |
4ae1c9f0 | 179 | fYieldArray[kGenPhi] = fPtPhi->Integral(fPtMin,fPtMax,(Double_t *)0,1.e-6); |
faf00b0d | 180 | #endif |
4ae1c9f0 | 181 | |
182 | // jpsi | |
183 | if(fHeaviestParticle<kGenJpsi)return; | |
184 | AliGenParam *genjpsi=0; | |
185 | Char_t nameJpsi[10]; | |
186 | snprintf(nameJpsi,10,"Jpsi"); | |
187 | genjpsi = new AliGenParam(fNPart, new AliGenEMlib(), AliGenEMlib::kJpsi, "DUMMY"); | |
188 | AddSource2Generator(nameJpsi,genjpsi); | |
189 | TF1 *fPtJpsi = genjpsi->GetPt(); | |
190 | #if ROOT_VERSION_CODE >= ROOT_VERSION(5,99,0) | |
191 | fYieldArray[kGenJpsi] = fPtJpsi->Integral(fPtMin,fPtMax,1.e-6); | |
192 | #else | |
193 | fYieldArray[kGenJpsi] = fPtJpsi->Integral(fPtMin,fPtMax,(Double_t *)0,1.e-6); | |
194 | #endif | |
195 | ||
e40b9538 | 196 | } |
197 | ||
198 | //------------------------------------------------------------------- | |
199 | void AliGenEMCocktail::AddSource2Generator(Char_t* nameSource, | |
200 | AliGenParam* const genSource) | |
201 | { | |
202 | // add sources to the cocktail | |
203 | Double_t phiMin = fPhiMin*180./TMath::Pi(); | |
204 | Double_t phiMax = fPhiMax*180./TMath::Pi(); | |
205 | ||
206 | genSource->SetPtRange(fPtMin, fPtMax); | |
207 | genSource->SetYRange(fYMin, fYMax); | |
208 | genSource->SetPhiRange(phiMin, phiMax); | |
209 | genSource->SetWeighting(fWeightingMode); | |
4ae1c9f0 | 210 | genSource->SetForceGammaConversion(fForceConv); |
e40b9538 | 211 | if (!gMC) genSource->SetDecayer(fDecayer); |
212 | genSource->Init(); | |
213 | ||
214 | AddGenerator(genSource,nameSource,1.); // Adding Generator | |
215 | } | |
216 | ||
217 | //------------------------------------------------------------------- | |
218 | void AliGenEMCocktail::Init() | |
219 | { | |
4ae1c9f0 | 220 | // Initialisation |
e40b9538 | 221 | TIter next(fEntries); |
222 | AliGenCocktailEntry *entry; | |
223 | if (fStack) { | |
224 | while((entry = (AliGenCocktailEntry*)next())) { | |
225 | entry->Generator()->SetStack(fStack); | |
226 | } | |
227 | } | |
228 | } | |
229 | ||
230 | //_________________________________________________________________________ | |
231 | void AliGenEMCocktail::Generate() | |
232 | { | |
4ae1c9f0 | 233 | // Generate event |
e40b9538 | 234 | TIter next(fEntries); |
235 | AliGenCocktailEntry *entry = 0; | |
236 | AliGenCocktailEntry *preventry = 0; | |
237 | AliGenerator* gen = 0; | |
238 | ||
239 | if (fHeader) delete fHeader; | |
240 | fHeader = new AliGenCocktailEventHeader("Electromagnetic Cocktail Header"); | |
241 | ||
242 | const TObjArray *partArray = gAlice->GetMCApp()->Particles(); | |
243 | ||
4ae1c9f0 | 244 | // Generate the vertex position used by all generators |
e40b9538 | 245 | if(fVertexSmear == kPerEvent) Vertex(); |
246 | ||
4ae1c9f0 | 247 | //Reseting stack |
e40b9538 | 248 | AliRunLoader * runloader = AliRunLoader::Instance(); |
249 | if (runloader) | |
250 | if (runloader->Stack()) | |
251 | runloader->Stack()->Clean(); | |
252 | ||
4ae1c9f0 | 253 | // Loop over generators and generate events |
e40b9538 | 254 | Int_t igen = 0; |
6078e216 | 255 | Float_t evPlane; |
256 | Rndm(&evPlane,1); | |
257 | evPlane*=TMath::Pi()*2; | |
e40b9538 | 258 | while((entry = (AliGenCocktailEntry*)next())) { |
259 | gen = entry->Generator(); | |
e40b9538 | 260 | gen->SetVertex(fVertex.At(0), fVertex.At(1), fVertex.At(2)); |
261 | ||
262 | if (fNPart > 0) { | |
263 | igen++; | |
264 | if (igen == 1) entry->SetFirst(0); | |
265 | else entry->SetFirst((partArray->GetEntriesFast())+1); | |
6078e216 | 266 | gen->SetEventPlane(evPlane); |
e40b9538 | 267 | gen->SetNumberParticles(fNPart); |
268 | gen->Generate(); | |
269 | entry->SetLast(partArray->GetEntriesFast()); | |
270 | preventry = entry; | |
271 | } | |
272 | } | |
273 | next.Reset(); | |
274 | ||
4ae1c9f0 | 275 | // Setting weights for proper absolute normalization |
e40b9538 | 276 | Int_t iPart, iMother; |
277 | Int_t pdgMother = 0; | |
278 | Double_t weight = 0.; | |
279 | Double_t dNdy = 0.; | |
280 | Int_t maxPart = partArray->GetEntriesFast(); | |
281 | for(iPart=0; iPart<maxPart; iPart++){ | |
282 | TParticle *part = gAlice->GetMCApp()->Particle(iPart); | |
283 | iMother = part->GetFirstMother(); | |
284 | TParticle *mother = 0; | |
285 | if (iMother>=0){ | |
286 | mother = gAlice->GetMCApp()->Particle(iMother); | |
287 | pdgMother = mother->GetPdgCode(); | |
288 | } | |
289 | else | |
290 | pdgMother = part->GetPdgCode(); | |
291 | switch (pdgMother){ | |
292 | case 111: | |
293 | dNdy = fYieldArray[kGenPizero]; | |
294 | break; | |
295 | case 221: | |
296 | dNdy = fYieldArray[kGenEta]; | |
297 | break; | |
298 | case 113: | |
299 | dNdy = fYieldArray[kGenRho]; | |
300 | break; | |
301 | case 223: | |
302 | dNdy = fYieldArray[kGenOmega]; | |
303 | break; | |
304 | case 331: | |
305 | dNdy = fYieldArray[kGenEtaprime]; | |
306 | break; | |
307 | case 333: | |
308 | dNdy = fYieldArray[kGenPhi]; | |
309 | break; | |
4ae1c9f0 | 310 | case 443: |
311 | dNdy = fYieldArray[kGenJpsi]; | |
312 | break; | |
e40b9538 | 313 | |
314 | default: | |
315 | dNdy = 0.; | |
316 | } | |
e40b9538 | 317 | weight = dNdy*part->GetWeight(); |
318 | part->SetWeight(weight); | |
319 | } | |
320 | ||
321 | fHeader->SetNProduced(maxPart); | |
322 | ||
323 | ||
324 | TArrayF eventVertex; | |
325 | eventVertex.Set(3); | |
326 | for (Int_t j=0; j < 3; j++) eventVertex[j] = fVertex[j]; | |
327 | ||
328 | fHeader->SetPrimaryVertex(eventVertex); | |
329 | ||
330 | gAlice->SetGenEventHeader(fHeader); | |
331 | } |