]>
Commit | Line | Data |
---|---|---|
415d9f5c | 1 | /************************************************************************** |
2 | * Copyright(c) 1998-2007, 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 | // Class for Kinematic Events | |
20 | // Author: Andreas Morsch, CERN | |
21 | //------------------------------------------------------------------------- | |
22 | #include <TArrow.h> | |
23 | #include <TMarker.h> | |
24 | #include <TH2F.h> | |
25 | #include <TTree.h> | |
26 | #include <TFile.h> | |
27 | #include <TParticle.h> | |
28 | #include <TClonesArray.h> | |
93df0e9b | 29 | #include <TRefArray.h> |
93836e1b | 30 | #include <TList.h> |
415d9f5c | 31 | |
32 | #include "AliLog.h" | |
33 | #include "AliMCEvent.h" | |
34 | #include "AliStack.h" | |
35 | #include "AliTrackReference.h" | |
36 | #include "AliHeader.h" | |
83626b6b | 37 | #include "AliGenEventHeader.h" |
415d9f5c | 38 | |
39 | ||
93836e1b | 40 | Int_t AliMCEvent::fgkBgLabelOffset(10000000); |
41 | ||
42 | ||
415d9f5c | 43 | AliMCEvent::AliMCEvent(): |
44 | AliVEvent(), | |
45 | fStack(0), | |
7aad0c47 | 46 | fMCParticles(0), |
415d9f5c | 47 | fMCParticleMap(0), |
fb9aeda5 | 48 | fHeader(new AliHeader()), |
93df0e9b | 49 | fTRBuffer(0), |
50 | fTrackReferences(new TClonesArray("AliTrackReference", 1000)), | |
415d9f5c | 51 | fTreeTR(0), |
52 | fTmpTreeTR(0), | |
53 | fTmpFileTR(0), | |
54 | fNprimaries(-1), | |
93836e1b | 55 | fNparticles(-1), |
56 | fSubsidiaryEvents(0), | |
57 | fPrimaryOffset(0), | |
7aad0c47 | 58 | fSecondaryOffset(0), |
59 | fExternal(0) | |
415d9f5c | 60 | { |
61 | // Default constructor | |
62 | } | |
63 | ||
64 | AliMCEvent::AliMCEvent(const AliMCEvent& mcEvnt) : | |
65 | AliVEvent(mcEvnt), | |
93836e1b | 66 | fStack(mcEvnt.fStack), |
67 | fMCParticles(mcEvnt.fMCParticles), | |
68 | fMCParticleMap(mcEvnt.fMCParticleMap), | |
69 | fHeader(mcEvnt.fHeader), | |
70 | fTRBuffer(mcEvnt.fTRBuffer), | |
71 | fTrackReferences(mcEvnt.fTrackReferences), | |
72 | fTreeTR(mcEvnt.fTreeTR), | |
73 | fTmpTreeTR(mcEvnt.fTmpTreeTR), | |
74 | fTmpFileTR(mcEvnt.fTmpFileTR), | |
75 | fNprimaries(mcEvnt.fNprimaries), | |
76 | fNparticles(mcEvnt.fNparticles), | |
77 | fSubsidiaryEvents(0), | |
78 | fPrimaryOffset(0), | |
7aad0c47 | 79 | fSecondaryOffset(0), |
80 | fExternal(0) | |
415d9f5c | 81 | { |
82 | // Copy constructor | |
83 | } | |
84 | ||
85 | ||
86 | AliMCEvent& AliMCEvent::operator=(const AliMCEvent& mcEvnt) | |
87 | { | |
88 | // assignment operator | |
89 | if (this!=&mcEvnt) { | |
90 | AliVEvent::operator=(mcEvnt); | |
91 | } | |
92 | ||
93 | return *this; | |
94 | } | |
95 | ||
96 | void AliMCEvent::ConnectTreeE (TTree* tree) | |
97 | { | |
98 | // Connect the event header tree | |
99 | tree->SetBranchAddress("Header", &fHeader); | |
100 | } | |
101 | ||
102 | void AliMCEvent::ConnectTreeK (TTree* tree) | |
103 | { | |
104 | // Connect the kinematics tree to the stack | |
7aad0c47 | 105 | if (!fMCParticles) fMCParticles = new TClonesArray("AliMCParticle",1000); |
106 | // | |
415d9f5c | 107 | fStack = fHeader->Stack(); |
108 | fStack->ConnectTree(tree); | |
109 | // | |
110 | // Load the event | |
111 | fStack->GetEvent(); | |
112 | fNparticles = fStack->GetNtrack(); | |
113 | fNprimaries = fStack->GetNprimary(); | |
ed97dc98 | 114 | Int_t iev = fHeader->GetEvent(); |
115 | Int_t ievr = fHeader->GetEventNrInRun(); | |
116 | ||
117 | AliInfo(Form("AliMCEvent# %5d %5d: Number of particles: %5d (all) %5d (primaries)\n", | |
118 | iev, ievr, fNparticles, fNprimaries)); | |
415d9f5c | 119 | |
120 | // This is a cache for the TParticles converted to MCParticles on user request | |
121 | if (fMCParticleMap) { | |
678a0069 | 122 | fMCParticleMap->Clear(); |
ee2774b7 | 123 | fMCParticles->Delete(); |
678a0069 | 124 | if (fNparticles>0) fMCParticleMap->Expand(fNparticles); |
125 | } | |
415d9f5c | 126 | else |
93df0e9b | 127 | fMCParticleMap = new TRefArray(fNparticles); |
415d9f5c | 128 | } |
129 | ||
130 | void AliMCEvent::ConnectTreeTR (TTree* tree) | |
131 | { | |
132 | // Connect the track reference tree | |
133 | fTreeTR = tree; | |
134 | ||
135 | if (fTreeTR->GetBranch("AliRun")) { | |
136 | if (fTmpFileTR) { | |
137 | fTmpFileTR->Close(); | |
138 | delete fTmpFileTR; | |
139 | } | |
140 | // This is an old format with one branch per detector not in synch with TreeK | |
141 | ReorderAndExpandTreeTR(); | |
142 | } else { | |
143 | // New format | |
93df0e9b | 144 | fTreeTR->SetBranchAddress("TrackReferences", &fTRBuffer); |
415d9f5c | 145 | } |
146 | } | |
147 | ||
148 | Int_t AliMCEvent::GetParticleAndTR(Int_t i, TParticle*& particle, TClonesArray*& trefs) | |
149 | { | |
150 | // Retrieve entry i | |
151 | if (i < 0 || i >= fNparticles) { | |
152 | AliWarning(Form("AliMCEventHandler::GetEntry: Index out of range")); | |
153 | particle = 0; | |
154 | trefs = 0; | |
155 | return (-1); | |
156 | } | |
157 | particle = fStack->Particle(i); | |
158 | if (fTreeTR) { | |
159 | fTreeTR->GetEntry(fStack->TreeKEntry(i)); | |
93df0e9b | 160 | trefs = fTRBuffer; |
415d9f5c | 161 | return trefs->GetEntries(); |
162 | } else { | |
163 | trefs = 0; | |
164 | return -1; | |
165 | } | |
166 | } | |
167 | ||
168 | ||
169 | void AliMCEvent::Clean() | |
170 | { | |
171 | // Clean-up before new trees are connected | |
fb9aeda5 | 172 | // if (fHeader) { |
173 | // delete fHeader; | |
174 | // fHeader = 0; | |
175 | // } | |
c6360f33 | 176 | |
177 | delete fStack; fStack = 0; | |
415d9f5c | 178 | |
179 | // Clear TR | |
93df0e9b | 180 | if (fTRBuffer) { |
ee2774b7 | 181 | fTRBuffer->Delete(); |
93df0e9b | 182 | delete fTRBuffer; |
183 | fTRBuffer = 0; | |
415d9f5c | 184 | } |
185 | } | |
186 | ||
187 | void AliMCEvent::FinishEvent() | |
188 | { | |
f1ac8a97 | 189 | // Clean-up after event |
c6360f33 | 190 | // |
191 | fStack->Reset(0); | |
192 | fMCParticles->Delete(); | |
678a0069 | 193 | fMCParticleMap->Clear(); |
568d5aed | 194 | if (fTRBuffer) |
195 | fTRBuffer->Delete(); | |
c6360f33 | 196 | fTrackReferences->Delete(); |
2bb794ba | 197 | fNparticles = -1; |
198 | fNprimaries = -1; | |
199 | fStack = 0; | |
93836e1b | 200 | // fSubsidiaryEvents->Clear(); |
201 | fSubsidiaryEvents = 0; | |
415d9f5c | 202 | } |
203 | ||
204 | ||
93df0e9b | 205 | |
415d9f5c | 206 | void AliMCEvent::DrawCheck(Int_t i, Int_t search) |
207 | { | |
208 | // | |
209 | // Simple event display for debugging | |
210 | if (!fTreeTR) { | |
211 | AliWarning("No Track Reference information available"); | |
212 | return; | |
213 | } | |
214 | ||
215 | if (i > -1 && i < fNparticles) { | |
216 | fTreeTR->GetEntry(fStack->TreeKEntry(i)); | |
217 | } else { | |
218 | AliWarning("AliMCEvent::GetEntry: Index out of range"); | |
219 | } | |
220 | ||
93df0e9b | 221 | Int_t nh = fTRBuffer->GetEntries(); |
415d9f5c | 222 | |
223 | ||
224 | if (search) { | |
225 | while(nh <= search && i < fNparticles - 1) { | |
226 | i++; | |
227 | fTreeTR->GetEntry(fStack->TreeKEntry(i)); | |
93df0e9b | 228 | nh = fTRBuffer->GetEntries(); |
415d9f5c | 229 | } |
230 | printf("Found Hits at %5d\n", i); | |
231 | } | |
232 | TParticle* particle = fStack->Particle(i); | |
233 | ||
234 | TH2F* h = new TH2F("", "", 100, -500, 500, 100, -500, 500); | |
235 | Float_t x0 = particle->Vx(); | |
236 | Float_t y0 = particle->Vy(); | |
237 | ||
238 | Float_t x1 = particle->Vx() + particle->Px() * 50.; | |
239 | Float_t y1 = particle->Vy() + particle->Py() * 50.; | |
240 | ||
241 | TArrow* a = new TArrow(x0, y0, x1, y1, 0.01); | |
242 | h->Draw(); | |
243 | a->SetLineColor(2); | |
244 | ||
245 | a->Draw(); | |
246 | ||
247 | for (Int_t ih = 0; ih < nh; ih++) { | |
93df0e9b | 248 | AliTrackReference* ref = (AliTrackReference*) fTRBuffer->At(ih); |
415d9f5c | 249 | TMarker* m = new TMarker(ref->X(), ref->Y(), 20); |
250 | m->Draw(); | |
251 | m->SetMarkerSize(0.4); | |
252 | ||
253 | } | |
254 | } | |
255 | ||
256 | ||
257 | void AliMCEvent::ReorderAndExpandTreeTR() | |
258 | { | |
259 | // | |
260 | // Reorder and expand the track reference tree in order to match the kinematics tree. | |
261 | // Copy the information from different branches into one | |
262 | // | |
263 | // TreeTR | |
264 | ||
265 | fTmpFileTR = new TFile("TrackRefsTmp.root", "recreate"); | |
266 | fTmpTreeTR = new TTree("TreeTR", "TrackReferences"); | |
93df0e9b | 267 | if (!fTRBuffer) fTRBuffer = new TClonesArray("AliTrackReference", 100); |
268 | fTmpTreeTR->Branch("TrackReferences", "TClonesArray", &fTRBuffer, 32000, 0); | |
415d9f5c | 269 | |
270 | ||
271 | // | |
272 | // Activate the used branches only. Otherwisw we get a bad memory leak. | |
273 | fTreeTR->SetBranchStatus("*", 0); | |
274 | fTreeTR->SetBranchStatus("AliRun.*", 1); | |
275 | fTreeTR->SetBranchStatus("ITS.*", 1); | |
276 | fTreeTR->SetBranchStatus("TPC.*", 1); | |
277 | fTreeTR->SetBranchStatus("TRD.*", 1); | |
278 | fTreeTR->SetBranchStatus("TOF.*", 1); | |
279 | fTreeTR->SetBranchStatus("FRAME.*", 1); | |
280 | fTreeTR->SetBranchStatus("MUON.*", 1); | |
281 | // | |
282 | // Connect the active branches | |
283 | TClonesArray* trefs[7]; | |
284 | for (Int_t i = 0; i < 7; i++) trefs[i] = 0; | |
285 | if (fTreeTR){ | |
286 | // make branch for central track references | |
287 | if (fTreeTR->GetBranch("AliRun")) fTreeTR->SetBranchAddress("AliRun", &trefs[0]); | |
288 | if (fTreeTR->GetBranch("ITS")) fTreeTR->SetBranchAddress("ITS", &trefs[1]); | |
289 | if (fTreeTR->GetBranch("TPC")) fTreeTR->SetBranchAddress("TPC", &trefs[2]); | |
290 | if (fTreeTR->GetBranch("TRD")) fTreeTR->SetBranchAddress("TRD", &trefs[3]); | |
291 | if (fTreeTR->GetBranch("TOF")) fTreeTR->SetBranchAddress("TOF", &trefs[4]); | |
292 | if (fTreeTR->GetBranch("FRAME")) fTreeTR->SetBranchAddress("FRAME", &trefs[5]); | |
293 | if (fTreeTR->GetBranch("MUON")) fTreeTR->SetBranchAddress("MUON", &trefs[6]); | |
294 | } | |
295 | ||
296 | Int_t np = fStack->GetNprimary(); | |
297 | Int_t nt = fTreeTR->GetEntries(); | |
298 | ||
299 | // | |
300 | // Loop over tracks and find the secondaries with the help of the kine tree | |
301 | Int_t ifills = 0; | |
302 | Int_t it = 0; | |
303 | Int_t itlast = 0; | |
304 | TParticle* part; | |
305 | ||
306 | for (Int_t ip = np - 1; ip > -1; ip--) { | |
307 | part = fStack->Particle(ip); | |
308 | // printf("Particle %5d %5d %5d %5d %5d %5d \n", | |
309 | // ip, part->GetPdgCode(), part->GetFirstMother(), part->GetFirstDaughter(), | |
310 | // part->GetLastDaughter(), part->TestBit(kTransportBit)); | |
311 | ||
312 | // Determine range of secondaries produced by this primary during transport | |
313 | Int_t dau1 = part->GetFirstDaughter(); | |
314 | if (dau1 < np) continue; // This particle has no secondaries produced during transport | |
315 | Int_t dau2 = -1; | |
316 | if (dau1 > -1) { | |
317 | Int_t inext = ip - 1; | |
318 | while (dau2 < 0) { | |
319 | if (inext >= 0) { | |
320 | part = fStack->Particle(inext); | |
321 | dau2 = part->GetFirstDaughter(); | |
322 | if (dau2 == -1 || dau2 < np) { | |
323 | dau2 = -1; | |
324 | } else { | |
325 | dau2--; | |
326 | } | |
327 | } else { | |
328 | dau2 = fStack->GetNtrack() - 1; | |
329 | } | |
330 | inext--; | |
331 | } // find upper bound | |
332 | } // dau2 < 0 | |
333 | ||
334 | ||
335 | // printf("Check (1) %5d %5d %5d %5d %5d \n", ip, np, it, dau1, dau2); | |
336 | // | |
337 | // Loop over reference hits and find secondary label | |
338 | // First the tricky part: find the entry in treeTR than contains the hits or | |
339 | // make sure that no hits exist. | |
340 | // | |
341 | Bool_t hasHits = kFALSE; | |
342 | Bool_t isOutside = kFALSE; | |
343 | ||
344 | it = itlast; | |
345 | while (!hasHits && !isOutside && it < nt) { | |
346 | fTreeTR->GetEntry(it++); | |
347 | for (Int_t ib = 0; ib < 7; ib++) { | |
348 | if (!trefs[ib]) continue; | |
349 | Int_t nh = trefs[ib]->GetEntries(); | |
350 | for (Int_t ih = 0; ih < nh; ih++) { | |
351 | AliTrackReference* tr = (AliTrackReference*) trefs[ib]->At(ih); | |
352 | Int_t label = tr->Label(); | |
353 | if (label >= dau1 && label <= dau2) { | |
354 | hasHits = kTRUE; | |
355 | itlast = it - 1; | |
356 | break; | |
357 | } | |
358 | if (label > dau2 || label < ip) { | |
359 | isOutside = kTRUE; | |
360 | itlast = it - 1; | |
361 | break; | |
362 | } | |
363 | } // hits | |
364 | if (hasHits || isOutside) break; | |
365 | } // branches | |
366 | } // entries | |
367 | ||
368 | if (!hasHits) { | |
369 | // Write empty entries | |
370 | for (Int_t id = dau1; (id <= dau2); id++) { | |
371 | fTmpTreeTR->Fill(); | |
372 | ifills++; | |
373 | } | |
374 | } else { | |
375 | // Collect all hits | |
376 | fTreeTR->GetEntry(itlast); | |
377 | for (Int_t id = dau1; (id <= dau2) && (dau1 > -1); id++) { | |
378 | for (Int_t ib = 0; ib < 7; ib++) { | |
379 | if (!trefs[ib]) continue; | |
380 | Int_t nh = trefs[ib]->GetEntries(); | |
381 | for (Int_t ih = 0; ih < nh; ih++) { | |
382 | AliTrackReference* tr = (AliTrackReference*) trefs[ib]->At(ih); | |
383 | Int_t label = tr->Label(); | |
384 | // Skip primaries | |
385 | if (label == ip) continue; | |
386 | if (label > dau2 || label < dau1) | |
387 | printf("AliMCEventHandler::Track Reference Label out of range !: %5d %5d %5d %5d \n", | |
388 | itlast, label, dau1, dau2); | |
389 | if (label == id) { | |
390 | // secondary found | |
391 | tr->SetDetectorId(ib-1); | |
93df0e9b | 392 | Int_t nref = fTRBuffer->GetEntriesFast(); |
393 | TClonesArray &lref = *fTRBuffer; | |
415d9f5c | 394 | new(lref[nref]) AliTrackReference(*tr); |
395 | } | |
396 | } // hits | |
397 | } // branches | |
398 | fTmpTreeTR->Fill(); | |
ee2774b7 | 399 | fTRBuffer->Delete(); |
415d9f5c | 400 | ifills++; |
401 | } // daughters | |
402 | } // has hits | |
403 | } // tracks | |
404 | ||
405 | // | |
406 | // Now loop again and write the primaries | |
407 | // | |
408 | it = nt - 1; | |
409 | for (Int_t ip = 0; ip < np; ip++) { | |
410 | Int_t labmax = -1; | |
411 | while (labmax < ip && it > -1) { | |
412 | fTreeTR->GetEntry(it--); | |
413 | for (Int_t ib = 0; ib < 7; ib++) { | |
414 | if (!trefs[ib]) continue; | |
415 | Int_t nh = trefs[ib]->GetEntries(); | |
416 | // | |
417 | // Loop over reference hits and find primary labels | |
418 | for (Int_t ih = 0; ih < nh; ih++) { | |
419 | AliTrackReference* tr = (AliTrackReference*) trefs[ib]->At(ih); | |
420 | Int_t label = tr->Label(); | |
421 | if (label < np && label > labmax) { | |
422 | labmax = label; | |
423 | } | |
424 | ||
425 | if (label == ip) { | |
426 | tr->SetDetectorId(ib-1); | |
93df0e9b | 427 | Int_t nref = fTRBuffer->GetEntriesFast(); |
428 | TClonesArray &lref = *fTRBuffer; | |
415d9f5c | 429 | new(lref[nref]) AliTrackReference(*tr); |
430 | } | |
431 | } // hits | |
432 | } // branches | |
433 | } // entries | |
434 | it++; | |
435 | fTmpTreeTR->Fill(); | |
ee2774b7 | 436 | fTRBuffer->Delete(); |
415d9f5c | 437 | ifills++; |
438 | } // tracks | |
439 | // Check | |
440 | ||
441 | ||
442 | // Clean-up | |
443 | delete fTreeTR; fTreeTR = 0; | |
444 | ||
445 | for (Int_t ib = 0; ib < 7; ib++) { | |
446 | if (trefs[ib]) { | |
447 | trefs[ib]->Clear(); | |
448 | delete trefs[ib]; | |
449 | trefs[ib] = 0; | |
450 | } | |
451 | } | |
452 | ||
453 | if (ifills != fStack->GetNtrack()) | |
454 | printf("AliMCEvent:Number of entries in TreeTR (%5d) unequal to TreeK (%5d) \n", | |
455 | ifills, fStack->GetNtrack()); | |
456 | ||
457 | fTmpTreeTR->Write(); | |
458 | fTreeTR = fTmpTreeTR; | |
459 | } | |
460 | ||
7aad0c47 | 461 | AliVParticle* AliMCEvent::GetTrack(Int_t i) const |
415d9f5c | 462 | { |
463 | // Get MC Particle i | |
93836e1b | 464 | // |
7aad0c47 | 465 | |
466 | if (fExternal) { | |
467 | return ((AliVParticle*) (fMCParticles->At(i))); | |
468 | } | |
469 | ||
93836e1b | 470 | // |
471 | // Check first if this explicitely accesses the subsidiary event | |
7aad0c47 | 472 | |
93836e1b | 473 | if (i >= BgLabelOffset()) { |
474 | if (fSubsidiaryEvents) { | |
475 | AliMCEvent* bgEvent = (AliMCEvent*) (fSubsidiaryEvents->At(1)); | |
476 | return (bgEvent->GetTrack(i - BgLabelOffset())); | |
477 | } else { | |
478 | return 0; | |
479 | } | |
480 | } | |
481 | ||
482 | // | |
93df0e9b | 483 | AliMCParticle *mcParticle = 0; |
484 | TParticle *particle = 0; | |
485 | TClonesArray *trefs = 0; | |
486 | Int_t ntref = 0; | |
487 | TRefArray *rarray = 0; | |
93836e1b | 488 | |
489 | ||
490 | ||
415d9f5c | 491 | // Out of range check |
492 | if (i < 0 || i >= fNparticles) { | |
493 | AliWarning(Form("AliMCEvent::GetEntry: Index out of range")); | |
494 | mcParticle = 0; | |
495 | return (mcParticle); | |
496 | } | |
93836e1b | 497 | |
415d9f5c | 498 | |
93836e1b | 499 | if (fSubsidiaryEvents) { |
500 | AliMCEvent* mc; | |
501 | Int_t idx = FindIndexAndEvent(i, mc); | |
502 | return (mc->GetTrack(idx)); | |
503 | } | |
415d9f5c | 504 | |
415d9f5c | 505 | // |
93836e1b | 506 | // First check If the MC Particle has been already cached |
415d9f5c | 507 | if(!fMCParticleMap->At(i)) { |
93df0e9b | 508 | // Get particle from the stack |
509 | particle = fStack->Particle(i); | |
510 | // Get track references from Tree TR | |
511 | if (fTreeTR) { | |
512 | fTreeTR->GetEntry(fStack->TreeKEntry(i)); | |
513 | trefs = fTRBuffer; | |
514 | ntref = trefs->GetEntriesFast(); | |
515 | rarray = new TRefArray(ntref); | |
516 | Int_t nen = fTrackReferences->GetEntriesFast(); | |
517 | for (Int_t j = 0; j < ntref; j++) { | |
518 | // Save the track references in a TClonesArray | |
519 | AliTrackReference* ref = dynamic_cast<AliTrackReference*>((*fTRBuffer)[j]); | |
520 | // Save the pointer in a TRefArray | |
521 | new ((*fTrackReferences)[nen]) AliTrackReference(*ref); | |
522 | rarray->AddAt((*fTrackReferences)[nen], j); | |
523 | nen++; | |
524 | } // loop over track references for entry i | |
525 | } // if TreeTR available | |
415d9f5c | 526 | Int_t nentries = fMCParticles->GetEntriesFast(); |
eee13e8d | 527 | new ((*fMCParticles)[nentries]) AliMCParticle(particle, rarray, i); |
415d9f5c | 528 | mcParticle = dynamic_cast<AliMCParticle*>((*fMCParticles)[nentries]); |
529 | fMCParticleMap->AddAt(mcParticle, i); | |
93836e1b | 530 | |
531 | TParticle* part = mcParticle->Particle(); | |
532 | Int_t imo = part->GetFirstMother(); | |
533 | Int_t id1 = part->GetFirstDaughter(); | |
534 | Int_t id2 = part->GetLastDaughter(); | |
535 | if (fPrimaryOffset > 0 || fSecondaryOffset > 0) { | |
536 | // Remapping of the mother and daughter indices | |
537 | if (imo < fNprimaries) { | |
538 | mcParticle->SetMother(imo + fPrimaryOffset); | |
539 | } else { | |
540 | mcParticle->SetMother(imo + fSecondaryOffset - fNprimaries); | |
541 | } | |
542 | ||
543 | if (id1 < fNprimaries) { | |
544 | mcParticle->SetFirstDaughter(id1 + fPrimaryOffset); | |
545 | mcParticle->SetLastDaughter (id2 + fPrimaryOffset); | |
546 | } else { | |
547 | mcParticle->SetFirstDaughter(id1 + fSecondaryOffset - fNprimaries); | |
548 | mcParticle->SetLastDaughter (id2 + fSecondaryOffset - fNprimaries); | |
549 | } | |
550 | ||
551 | ||
552 | if (i > fNprimaries) { | |
553 | mcParticle->SetLabel(i + fPrimaryOffset); | |
554 | } else { | |
555 | mcParticle->SetLabel(i + fSecondaryOffset - fNprimaries); | |
556 | } | |
557 | ||
558 | ||
559 | } else { | |
560 | mcParticle->SetFirstDaughter(id1); | |
561 | mcParticle->SetLastDaughter (id2); | |
562 | mcParticle->SetMother (imo); | |
563 | } | |
564 | ||
415d9f5c | 565 | } else { |
566 | mcParticle = dynamic_cast<AliMCParticle*>(fMCParticleMap->At(i)); | |
567 | } | |
93836e1b | 568 | |
569 | ||
415d9f5c | 570 | return mcParticle; |
571 | } | |
572 | ||
93836e1b | 573 | AliGenEventHeader* AliMCEvent::GenEventHeader() {return (fHeader->GenEventHeader());} |
574 | ||
575 | ||
576 | void AliMCEvent::AddSubsidiaryEvent(AliMCEvent* event) | |
577 | { | |
578 | // Add a subsidiary event to the list; for example merged background event. | |
579 | if (!fSubsidiaryEvents) { | |
580 | TList* events = new TList(); | |
581 | events->Add(new AliMCEvent(*this)); | |
582 | fSubsidiaryEvents = events; | |
583 | } | |
584 | ||
585 | fSubsidiaryEvents->Add(event); | |
586 | } | |
587 | ||
588 | Int_t AliMCEvent::FindIndexAndEvent(Int_t oldidx, AliMCEvent*& event) const | |
589 | { | |
590 | // Find the index and event in case of composed events like signal + background | |
591 | TIter next(fSubsidiaryEvents); | |
592 | next.Reset(); | |
593 | if (oldidx < fNprimaries) { | |
594 | while((event = (AliMCEvent*)next())) { | |
595 | if (oldidx < (event->GetPrimaryOffset() + event->GetNumberOfPrimaries())) break; | |
596 | } | |
597 | return (oldidx - event->GetPrimaryOffset()); | |
598 | } else { | |
599 | while((event = (AliMCEvent*)next())) { | |
600 | if (oldidx < (event->GetSecondaryOffset() + (event->GetNumberOfTracks() - event->GetNumberOfPrimaries()))) break; | |
601 | } | |
602 | return (oldidx - event->GetSecondaryOffset() + event->GetNumberOfPrimaries()); | |
603 | } | |
604 | } | |
605 | ||
606 | Int_t AliMCEvent::BgLabelToIndex(Int_t label) | |
607 | { | |
608 | // Convert a background label to an absolute index | |
609 | if (fSubsidiaryEvents) { | |
610 | AliMCEvent* bgEvent = (AliMCEvent*) (fSubsidiaryEvents->At(1)); | |
611 | label -= BgLabelOffset(); | |
612 | if (label < bgEvent->GetNumberOfPrimaries()) { | |
613 | label += bgEvent->GetPrimaryOffset(); | |
614 | } else { | |
615 | label += (bgEvent->GetSecondaryOffset() - fNprimaries); | |
616 | } | |
617 | } | |
618 | return (label); | |
619 | } | |
620 | ||
621 | ||
622 | Bool_t AliMCEvent::IsPhysicalPrimary(Int_t i) | |
623 | { | |
624 | // | |
625 | // Delegate to subevent if necesarry | |
626 | ||
627 | ||
628 | if (!fSubsidiaryEvents) { | |
629 | return fStack->IsPhysicalPrimary(i); | |
630 | } else { | |
631 | AliMCEvent* evt = 0; | |
632 | Int_t idx = FindIndexAndEvent(i, evt); | |
633 | return (evt->IsPhysicalPrimary(idx)); | |
634 | } | |
635 | } | |
636 | ||
637 | void AliMCEvent::InitEvent() | |
638 | { | |
35152e4b | 639 | // |
640 | // Initialize the subsidiary event structure | |
93836e1b | 641 | if (fSubsidiaryEvents) { |
642 | TIter next(fSubsidiaryEvents); | |
643 | AliMCEvent* evt; | |
644 | fNprimaries = 0; | |
645 | fNparticles = 0; | |
646 | ||
647 | while((evt = (AliMCEvent*)next())) { | |
648 | fNprimaries += evt->GetNumberOfPrimaries(); | |
649 | fNparticles += evt->GetNumberOfTracks(); | |
650 | } | |
651 | ||
652 | Int_t ioffp = 0; | |
653 | Int_t ioffs = fNprimaries; | |
654 | next.Reset(); | |
655 | ||
656 | while((evt = (AliMCEvent*)next())) { | |
657 | evt->SetPrimaryOffset(ioffp); | |
658 | evt->SetSecondaryOffset(ioffs); | |
659 | ioffp += evt->GetNumberOfPrimaries(); | |
660 | ioffs += (evt->GetNumberOfTracks() - evt->GetNumberOfPrimaries()); | |
661 | } | |
662 | } | |
663 | } | |
35152e4b | 664 | void AliMCEvent::PreReadAll() |
665 | { | |
666 | // Preread the MC information | |
667 | Int_t i; | |
668 | // secondaries | |
669 | for (i = fStack->GetNprimary(); i < fStack->GetNtrack(); i++) | |
670 | { | |
671 | GetTrack(i); | |
672 | } | |
673 | // primaries | |
674 | for (i = 0; i < fStack->GetNprimary(); i++) | |
675 | { | |
676 | GetTrack(i); | |
677 | } | |
678 | ||
679 | ||
680 | } | |
93836e1b | 681 | |
415d9f5c | 682 | |
683 | ClassImp(AliMCEvent) |