PreReading of MC information on demand.
[u/mrichter/AliRoot.git] / STEER / AliMCEventHandler.cxx
CommitLineData
93836e1b 1/**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved *
5fe09262 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// Class AliMCEventHandler
19// This class gives access to MC truth during the analysis.
20// Monte Carlo truth is containe in the kinematics tree (produced particles) and
21// the tree of reference hits.
22//
23// Origin: Andreas Morsch, CERN, andreas.morsch@cern.ch
24//---------------------------------------------------------------------------------
25
26
27
28#include "AliMCEventHandler.h"
415d9f5c 29#include "AliMCEvent.h"
93836e1b 30#include "AliMCParticle.h"
82418625 31#include "AliPDG.h"
5fe09262 32#include "AliTrackReference.h"
33#include "AliHeader.h"
34#include "AliStack.h"
2c36081e 35#include "AliLog.h"
5fe09262 36
37#include <TTree.h>
38#include <TFile.h>
93836e1b 39#include <TList.h>
5fe09262 40#include <TParticle.h>
0a05cd41 41#include <TString.h>
5fe09262 42#include <TClonesArray.h>
43#include <TDirectoryFile.h>
5fe09262 44
45ClassImp(AliMCEventHandler)
46
47AliMCEventHandler::AliMCEventHandler() :
d2f1d9ef 48 AliVEventHandler(),
415d9f5c 49 fMCEvent(new AliMCEvent()),
5fe09262 50 fFileE(0),
51 fFileK(0),
52 fFileTR(0),
5fe09262 53 fTreeE(0),
54 fTreeK(0),
55 fTreeTR(0),
5efedd31 56 fDirK(0),
57 fDirTR(0),
da97a08a 58 fParticleSelected(0),
59 fLabelMap(0),
5fe09262 60 fNEvent(-1),
61 fEvent(-1),
0a05cd41 62 fPathName(new TString("./")),
9aea8469 63 fExtension(""),
64 fFileNumber(0),
969c7896 65 fEventsPerFile(0),
c8b7b5d3 66 fReadTR(kTRUE),
93836e1b 67 fInitOk(kFALSE),
68 fSubsidiaryHandlers(0),
35152e4b 69 fEventsInContainer(0),
70 fPreReadMode(kNoPreRead)
5fe09262 71{
82418625 72 //
73 // Default constructor
74 //
75 // Be sure to add all particles to the PDG database
76 AliPDG::AddParticlesToPdgDataBase();
5fe09262 77}
78
79AliMCEventHandler::AliMCEventHandler(const char* name, const char* title) :
d2f1d9ef 80 AliVEventHandler(name, title),
415d9f5c 81 fMCEvent(new AliMCEvent()),
5fe09262 82 fFileE(0),
83 fFileK(0),
84 fFileTR(0),
5fe09262 85 fTreeE(0),
86 fTreeK(0),
87 fTreeTR(0),
5efedd31 88 fDirK(0),
89 fDirTR(0),
da97a08a 90 fParticleSelected(0),
91 fLabelMap(0),
5fe09262 92 fNEvent(-1),
93 fEvent(-1),
0a05cd41 94 fPathName(new TString("./")),
9aea8469 95 fExtension(""),
96 fFileNumber(0),
0cd61c1d 97 fEventsPerFile(0),
c8b7b5d3 98 fReadTR(kTRUE),
93836e1b 99 fInitOk(kFALSE),
100 fSubsidiaryHandlers(0),
35152e4b 101 fEventsInContainer(0),
102 fPreReadMode(kNoPreRead)
5fe09262 103{
82418625 104 //
105 // Constructor
106 //
107 // Be sure to add all particles to the PDG database
108 AliPDG::AddParticlesToPdgDataBase();
5fe09262 109}
110AliMCEventHandler::~AliMCEventHandler()
111{
112 // Destructor
415d9f5c 113 delete fMCEvent;
5fe09262 114 delete fFileE;
115 delete fFileK;
116 delete fFileTR;
117}
118
300d5701 119Bool_t AliMCEventHandler::Init(Option_t* opt)
5fe09262 120{
121 // Initialize input
122 //
6073f8c9 123 if (!(strcmp(opt, "proof")) || !(strcmp(opt, "local"))) return kTRUE;
124 //
0a05cd41 125 fFileE = TFile::Open(Form("%sgalice.root", fPathName->Data()));
c8b7b5d3 126 if (!fFileE) {
127 AliError(Form("AliMCEventHandler:galice.root not found in directory %s ! \n", fPathName->Data()));
128 fInitOk = kFALSE;
129 return kFALSE;
130 }
131
415d9f5c 132 //
133 // Tree E
5fe09262 134 fFileE->GetObject("TE", fTreeE);
415d9f5c 135 // Connect Tree E to the MCEvent
136 fMCEvent->ConnectTreeE(fTreeE);
5fe09262 137 fNEvent = fTreeE->GetEntries();
138 //
139 // Tree K
0a05cd41 140 fFileK = TFile::Open(Form("%sKinematics%s.root", fPathName->Data(), fExtension));
c8b7b5d3 141 if (!fFileK) {
142 AliError(Form("AliMCEventHandler:Kinematics.root not found in directory %s ! \n", fPathName));
143 fInitOk = kFALSE;
b72029a3 144 return kTRUE;
c8b7b5d3 145 }
146
9aea8469 147 fEventsPerFile = fFileK->GetNkeys() - fFileK->GetNProcessIDs();
5fe09262 148 //
149 // Tree TR
969c7896 150 if (fReadTR) {
c8b7b5d3 151 fFileTR = TFile::Open(Form("%sTrackRefs%s.root", fPathName->Data(), fExtension));
152 if (!fFileTR) {
153 AliError(Form("AliMCEventHandler:TrackRefs.root not found in directory %s ! \n", fPathName->Data()));
154 fInitOk = kFALSE;
b72029a3 155 return kTRUE;
c8b7b5d3 156 }
969c7896 157 }
5fe09262 158 //
159 // Reset the event number
2c36081e 160 fEvent = -1;
161 fFileNumber = 0;
b544c64d 162 AliInfo(Form("Number of events in this directory %5d \n", fNEvent));
c8b7b5d3 163 fInitOk = kTRUE;
93836e1b 164
165
166 if (fSubsidiaryHandlers) {
167 TIter next(fSubsidiaryHandlers);
168 AliMCEventHandler *handler;
169 while((handler = (AliMCEventHandler*)next())) {
170 handler->Init(opt);
171 handler->SetNumberOfEventsInContainer(fNEvent);
172 }
173 }
174
5fe09262 175 return kTRUE;
5fe09262 176}
177
9aea8469 178Bool_t AliMCEventHandler::GetEvent(Int_t iev)
179{
180 // Load the event number iev
2c36081e 181 //
182 // Calculate the file number
c8b7b5d3 183 if (!fInitOk) return kFALSE;
184
415d9f5c 185 Int_t inew = iev / fEventsPerFile;
9aea8469 186 if (inew != fFileNumber) {
187 fFileNumber = inew;
188 if (!OpenFile(fFileNumber)){
189 return kFALSE;
190 }
191 }
2c36081e 192 // Folder name
5fe09262 193 char folder[20];
9aea8469 194 sprintf(folder, "Event%d", iev);
5fe09262 195 // TreeE
9aea8469 196 fTreeE->GetEntry(iev);
5fe09262 197 // Tree K
5efedd31 198 fFileK->GetObject(folder, fDirK);
199 if (!fDirK) {
2c36081e 200 AliWarning(Form("AliMCEventHandler: Event #%5d not found\n", iev));
9aea8469 201 return kFALSE;
202 }
b72029a3 203
5efedd31 204 fDirK ->GetObject("TreeK", fTreeK);
415d9f5c 205 // Connect TreeK to MCEvent
206 fMCEvent->ConnectTreeK(fTreeK);
5fe09262 207 //Tree TR
07bd4750 208 if (fFileTR) {
5efedd31 209 // Check which format has been read
210 fFileTR->GetObject(folder, fDirTR);
211 fDirTR->GetObject("TreeTR", fTreeTR);
415d9f5c 212 //
213 // Connect TR to MCEvent
214 fMCEvent->ConnectTreeTR(fTreeTR);
5fe09262 215 }
93836e1b 216
5fe09262 217 //
5fe09262 218 return kTRUE;
9aea8469 219}
220
221Bool_t AliMCEventHandler::OpenFile(Int_t i)
222{
223 // Open file i
9aea8469 224 if (i > 0) {
225 fExtension = Form("%d", i);
226 } else {
227 fExtension = "";
228 }
229
230
231 delete fFileK;
0a05cd41 232 fFileK = TFile::Open(Form("%sKinematics%s.root", fPathName->Data(), fExtension));
9aea8469 233 if (!fFileK) {
c8b7b5d3 234 AliError(Form("AliMCEventHandler:Kinematics%s.root not found in directory %s ! \n", fExtension, fPathName->Data()));
235 fInitOk = kFALSE;
236 return kFALSE;
9aea8469 237 }
238
c8b7b5d3 239 if (fReadTR) {
240 delete fFileTR;
241 fFileTR = TFile::Open(Form("%sTrackRefs%s.root", fPathName->Data(), fExtension));
242 if (!fFileTR) {
243 AliWarning(Form("AliMCEventHandler:TrackRefs%s.root not found in directory %s ! \n", fExtension, fPathName->Data()));
244 fInitOk = kFALSE;
245 return kFALSE;
246 }
9aea8469 247 }
248
c8b7b5d3 249 fInitOk = kTRUE;
93836e1b 250
c8b7b5d3 251 return kTRUE;
9aea8469 252}
253
ed97dc98 254Bool_t AliMCEventHandler::BeginEvent(Long64_t entry)
9aea8469 255{
da97a08a 256 fParticleSelected.Delete();
257 fLabelMap.Delete();
9aea8469 258 // Read the next event
93836e1b 259
260 if (fEventsInContainer != 0) {
11e4d716 261 entry = (Long64_t) ( entry * Float_t(fNEvent) / Float_t (fEventsInContainer));
93836e1b 262 }
263
264
ed97dc98 265 if (entry == -1) {
266 fEvent++;
267 entry = fEvent;
268 } else {
269 fEvent = entry;
270 }
271
272 if (entry >= fNEvent) {
93836e1b 273 AliWarning(Form("AliMCEventHandler: Event number out of range %5d %5d\n", entry, fNEvent));
9aea8469 274 return kFALSE;
275 }
93836e1b 276
277 Bool_t result = GetEvent(entry);
278
279 if (fSubsidiaryHandlers) {
280 TIter next(fSubsidiaryHandlers);
281 AliMCEventHandler *handler;
282 while((handler = (AliMCEventHandler*)next())) {
283 handler->BeginEvent(entry);
284 }
285 next.Reset();
286 while((handler = (AliMCEventHandler*)next())) {
287 fMCEvent->AddSubsidiaryEvent(handler->MCEvent());
288 }
289 fMCEvent->InitEvent();
290 }
35152e4b 291
292 if (fPreReadMode == kLmPreRead) {
293 fMCEvent->PreReadAll();
294 }
295
93836e1b 296 return result;
297
5fe09262 298}
299
da97a08a 300void AliMCEventHandler::SelectParticle(Int_t i){
301 // taking the absolute values here, need to take care
302 // of negative daughter and mother
303 // IDs when setting!
93836e1b 304 if (TMath::Abs(i) >= AliMCEvent::BgLabelOffset()) i = fMCEvent->BgLabelToIndex(TMath::Abs(i));
305 if(!IsParticleSelected(TMath::Abs(i)))fParticleSelected.Add(TMath::Abs(i),1);
da97a08a 306}
307
308Bool_t AliMCEventHandler::IsParticleSelected(Int_t i) {
309 // taking the absolute values here, need to take
310 // care with negative daughter and mother
311 // IDs when setting!
312 return (fParticleSelected.GetValue(TMath::Abs(i))==1);
313}
314
315
316void AliMCEventHandler::CreateLabelMap(){
317
318 //
319 // this should be called once all selections where done
320 //
321
322 fLabelMap.Delete();
323 if(!fMCEvent){
324 fParticleSelected.Delete();
325 return;
326 }
327
328 VerifySelectedParticles();
da97a08a 329
330 Int_t iNew = 0;
93836e1b 331 for(int i = 0;i < fMCEvent->GetNumberOfTracks();++i){
da97a08a 332 if(IsParticleSelected(i)){
333 fLabelMap.Add(i,iNew);
334 iNew++;
335 }
336 }
337}
338
339Int_t AliMCEventHandler::GetNewLabel(Int_t i) {
93836e1b 340 // Gets the label from the new created Map
da97a08a 341 // Call CreatLabelMap before
342 // otherwise only 0 returned
343 return fLabelMap.GetValue(TMath::Abs(i));
344}
345
346void AliMCEventHandler::VerifySelectedParticles(){
347
348 //
349 // Make sure that each particle has at least it's predecessors
350 // selected so that we have the complete ancestry tree
351 // Private, should be only called by CreateLabelMap
352
353 if(!fMCEvent){
93836e1b 354 fParticleSelected.Delete();
355 return;
da97a08a 356 }
da97a08a 357
93836e1b 358 Int_t nprim = fMCEvent->GetNumberOfPrimaries();
359
360 for(int i = 0;i < fMCEvent->GetNumberOfTracks(); ++i){
361 if(i < nprim){
362 SelectParticle(i);// take all primaries
363 continue;
364 }
365
366 if(!IsParticleSelected(i))continue;
367
63a1afff 368 AliMCParticle* mcpart = (AliMCParticle*) fMCEvent->GetTrack(i);
93836e1b 369 Int_t imo = mcpart->GetMother();
370 while((imo >= nprim)&&!IsParticleSelected(imo)){
371 // Mother not yet selected
372 SelectParticle(imo);
0f0c06de 373 mcpart = (AliMCParticle*) fMCEvent->GetTrack(imo);
93836e1b 374 imo = mcpart->GetMother();
375 }
376 // after last step we may have an unselected primary
da97a08a 377 // mother
378 if(imo>=0){
379 if(!IsParticleSelected(imo))
380 SelectParticle(imo);
381 }
382 }// loop over all tracks
383}
384
5fe09262 385Int_t AliMCEventHandler::GetParticleAndTR(Int_t i, TParticle*& particle, TClonesArray*& trefs)
386{
387 // Retrieve entry i
c8b7b5d3 388 if (!fInitOk) {
389 return 0;
390 } else {
391 return (fMCEvent->GetParticleAndTR(i, particle, trefs));
392 }
5fe09262 393}
394
2d8f26f6 395void AliMCEventHandler::DrawCheck(Int_t i, Int_t search)
5fe09262 396{
397 // Retrieve entry i and draw momentum vector and hits
415d9f5c 398 fMCEvent->DrawCheck(i, search);
5fe09262 399}
400
890126ab 401Bool_t AliMCEventHandler::Notify(const char *path)
5fe09262 402{
53faeca4 403 // Notify about directory change
404 // The directory is taken from the 'path' argument
405 // Reconnect trees
406 TString fileName(path);
407 if(fileName.Contains("AliESDs.root")){
408 fileName.ReplaceAll("AliESDs.root", "");
409 }
c8b5ce3a 410 else if(fileName.Contains("AliAOD.root")){
411 fileName.ReplaceAll("AliAOD.root", "");
412 }
53faeca4 413 else if(fileName.Contains("galice.root")){
414 // for running with galice and kinematics alone...
415 fileName.ReplaceAll("galice.root", "");
416 }
39734d21 417 else if (fileName.BeginsWith("root:")) {
418 fileName.Append("?ZIP=");
419 }
93836e1b 420
53faeca4 421 *fPathName = fileName;
93836e1b 422 AliInfo(Form("Path: -%s-\n", fPathName->Data()));
53faeca4 423
f6065654 424 ResetIO();
425 InitIO("");
426
93836e1b 427// Handle subsidiary handlers
428 if (fSubsidiaryHandlers) {
429 TIter next(fSubsidiaryHandlers);
430 AliMCEventHandler *handler;
431 while((handler = (AliMCEventHandler*) next())) {
432 TString* spath = handler->GetInputPath();
433 if (spath->Contains("merged")) {
434 if (! fPathName->IsNull()) {
435 handler->Notify(Form("%s/../.", fPathName->Data()));
436 } else {
437 handler->Notify("../");
438 }
439 }
440 }
441 }
442
5fe09262 443 return kTRUE;
444}
36e82a52 445
5fe09262 446void AliMCEventHandler::ResetIO()
447{
5efedd31 448// Clear header and stack
b72029a3 449
aff49450 450 if (fInitOk) fMCEvent->Clean();
5efedd31 451
36e82a52 452// Delete Tree E
415d9f5c 453 delete fTreeE; fTreeE = 0;
b72029a3 454
5efedd31 455// Reset files
89f249fc 456 if (fFileE) {delete fFileE; fFileE = 0;}
457 if (fFileK) {delete fFileK; fFileK = 0;}
458 if (fFileTR) {delete fFileTR; fFileTR = 0;}
19cfde04 459 fExtension="";
c8b7b5d3 460 fInitOk = kFALSE;
93836e1b 461
462 if (fSubsidiaryHandlers) {
463 TIter next(fSubsidiaryHandlers);
464 AliMCEventHandler *handler;
465 while((handler = (AliMCEventHandler*)next())) {
466 handler->ResetIO();
467 }
468 }
469
5fe09262 470}
471
472
473Bool_t AliMCEventHandler::FinishEvent()
474{
5efedd31 475 // Clean-up after each event
476 delete fDirTR; fDirTR = 0;
477 delete fDirK; fDirK = 0;
aff49450 478 if (fInitOk) fMCEvent->FinishEvent();
93836e1b 479
480 if (fSubsidiaryHandlers) {
481 TIter next(fSubsidiaryHandlers);
482 AliMCEventHandler *handler;
483 while((handler = (AliMCEventHandler*)next())) {
484 handler->FinishEvent();
485 }
486 }
487
5fe09262 488 return kTRUE;
489}
490
491Bool_t AliMCEventHandler::Terminate()
492{
493 // Dummy
494 return kTRUE;
495}
496
497Bool_t AliMCEventHandler::TerminateIO()
498{
499 // Dummy
500 return kTRUE;
501}
502
0a05cd41 503
0931e76a 504void AliMCEventHandler::SetInputPath(const char* fname)
0a05cd41 505{
506 // Set the input path name
507 delete fPathName;
508 fPathName = new TString(fname);
509}
93836e1b 510
511void AliMCEventHandler::AddSubsidiaryHandler(AliMCEventHandler* handler)
512{
513 // Add a subsidiary handler. For example for background events
514
515 if (!fSubsidiaryHandlers) fSubsidiaryHandlers = new TList();
516 fSubsidiaryHandlers->Add(handler);
517}