Moving the classes that belong to the following libraries: STEERBase, ESD, CDB, AOD...
[u/mrichter/AliRoot.git] / STEER / STEERBase / 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
3f61d763 113 delete fPathName;
415d9f5c 114 delete fMCEvent;
5fe09262 115 delete fFileE;
116 delete fFileK;
117 delete fFileTR;
118}
119
300d5701 120Bool_t AliMCEventHandler::Init(Option_t* opt)
5fe09262 121{
122 // Initialize input
123 //
6073f8c9 124 if (!(strcmp(opt, "proof")) || !(strcmp(opt, "local"))) return kTRUE;
125 //
0a05cd41 126 fFileE = TFile::Open(Form("%sgalice.root", fPathName->Data()));
c8b7b5d3 127 if (!fFileE) {
128 AliError(Form("AliMCEventHandler:galice.root not found in directory %s ! \n", fPathName->Data()));
129 fInitOk = kFALSE;
130 return kFALSE;
131 }
132
415d9f5c 133 //
134 // Tree E
5fe09262 135 fFileE->GetObject("TE", fTreeE);
415d9f5c 136 // Connect Tree E to the MCEvent
137 fMCEvent->ConnectTreeE(fTreeE);
5fe09262 138 fNEvent = fTreeE->GetEntries();
139 //
140 // Tree K
0a05cd41 141 fFileK = TFile::Open(Form("%sKinematics%s.root", fPathName->Data(), fExtension));
c8b7b5d3 142 if (!fFileK) {
65b25288 143 AliError(Form("AliMCEventHandler:Kinematics.root not found in directory %s ! \n", fPathName->Data()));
c8b7b5d3 144 fInitOk = kFALSE;
b72029a3 145 return kTRUE;
c8b7b5d3 146 }
147
9aea8469 148 fEventsPerFile = fFileK->GetNkeys() - fFileK->GetNProcessIDs();
5fe09262 149 //
150 // Tree TR
969c7896 151 if (fReadTR) {
c8b7b5d3 152 fFileTR = TFile::Open(Form("%sTrackRefs%s.root", fPathName->Data(), fExtension));
153 if (!fFileTR) {
154 AliError(Form("AliMCEventHandler:TrackRefs.root not found in directory %s ! \n", fPathName->Data()));
155 fInitOk = kFALSE;
b72029a3 156 return kTRUE;
c8b7b5d3 157 }
969c7896 158 }
5fe09262 159 //
160 // Reset the event number
2c36081e 161 fEvent = -1;
162 fFileNumber = 0;
b544c64d 163 AliInfo(Form("Number of events in this directory %5d \n", fNEvent));
c8b7b5d3 164 fInitOk = kTRUE;
93836e1b 165
166
167 if (fSubsidiaryHandlers) {
168 TIter next(fSubsidiaryHandlers);
169 AliMCEventHandler *handler;
170 while((handler = (AliMCEventHandler*)next())) {
171 handler->Init(opt);
172 handler->SetNumberOfEventsInContainer(fNEvent);
173 }
174 }
175
5fe09262 176 return kTRUE;
5fe09262 177}
178
9aea8469 179Bool_t AliMCEventHandler::GetEvent(Int_t iev)
180{
181 // Load the event number iev
2c36081e 182 //
183 // Calculate the file number
68793909 184 if (!fInitOk) return kFALSE;
c8b7b5d3 185
68793909 186 Int_t inew = iev / fEventsPerFile;
187 if (inew != fFileNumber) {
188 fFileNumber = inew;
189 if (!OpenFile(fFileNumber)){
190 return kFALSE;
9aea8469 191 }
68793909 192 }
193 // Folder name
194 char folder[20];
76880c9a 195 snprintf(folder, 20, "Event%d", iev);
68793909 196 // TreeE
197 fTreeE->GetEntry(iev);
198 // Tree K
199 fFileK->GetObject(folder, fDirK);
200 if (!fDirK) {
201 AliWarning(Form("AliMCEventHandler: Event #%5d - Cannot get kinematics\n", iev));
202 return kFALSE;
203 }
b72029a3 204
68793909 205 fDirK ->GetObject("TreeK", fTreeK);
206 if (!fTreeK) {
207 AliError(Form("AliMCEventHandler: Event #%5d - Cannot get TreeK\n",iev));
208 return kFALSE;
209 }
210 // Connect TreeK to MCEvent
211 fMCEvent->ConnectTreeK(fTreeK);
212 //Tree TR
213 if (fFileTR) {
214 // Check which format has been read
215 fFileTR->GetObject(folder, fDirTR);
216 if (!fDirTR) {
217 AliError(Form("AliMCEventHandler: Event #%5d - Cannot get track references\n",iev));
218 return kFALSE;
219 }
220
221 fDirTR->GetObject("TreeTR", fTreeTR);
5fe09262 222 //
68793909 223 if (!fTreeTR) {
224 AliError(Form("AliMCEventHandler: Event #%5d - Cannot get TreeTR\n",iev));
225 return kFALSE;
226 }
227 // Connect TR to MCEvent
228 fMCEvent->ConnectTreeTR(fTreeTR);
229 }
230
231 //
232 return kTRUE;
9aea8469 233}
234
235Bool_t AliMCEventHandler::OpenFile(Int_t i)
236{
237 // Open file i
9aea8469 238 if (i > 0) {
239 fExtension = Form("%d", i);
240 } else {
241 fExtension = "";
242 }
243
244
245 delete fFileK;
0a05cd41 246 fFileK = TFile::Open(Form("%sKinematics%s.root", fPathName->Data(), fExtension));
9aea8469 247 if (!fFileK) {
c8b7b5d3 248 AliError(Form("AliMCEventHandler:Kinematics%s.root not found in directory %s ! \n", fExtension, fPathName->Data()));
249 fInitOk = kFALSE;
250 return kFALSE;
9aea8469 251 }
252
c8b7b5d3 253 if (fReadTR) {
254 delete fFileTR;
255 fFileTR = TFile::Open(Form("%sTrackRefs%s.root", fPathName->Data(), fExtension));
256 if (!fFileTR) {
257 AliWarning(Form("AliMCEventHandler:TrackRefs%s.root not found in directory %s ! \n", fExtension, fPathName->Data()));
258 fInitOk = kFALSE;
259 return kFALSE;
260 }
9aea8469 261 }
262
c8b7b5d3 263 fInitOk = kTRUE;
93836e1b 264
c8b7b5d3 265 return kTRUE;
9aea8469 266}
267
ed97dc98 268Bool_t AliMCEventHandler::BeginEvent(Long64_t entry)
9aea8469 269{
aa7e002c 270 // Begin event
da97a08a 271 fParticleSelected.Delete();
272 fLabelMap.Delete();
9aea8469 273 // Read the next event
93836e1b 274
275 if (fEventsInContainer != 0) {
11e4d716 276 entry = (Long64_t) ( entry * Float_t(fNEvent) / Float_t (fEventsInContainer));
93836e1b 277 }
278
279
ed97dc98 280 if (entry == -1) {
281 fEvent++;
282 entry = fEvent;
283 } else {
284 fEvent = entry;
285 }
286
287 if (entry >= fNEvent) {
73bbf779 288 AliWarning(Form("AliMCEventHandler: Event number out of range %5lld %5d\n", entry, fNEvent));
9aea8469 289 return kFALSE;
290 }
93836e1b 291
292 Bool_t result = GetEvent(entry);
293
294 if (fSubsidiaryHandlers) {
295 TIter next(fSubsidiaryHandlers);
296 AliMCEventHandler *handler;
297 while((handler = (AliMCEventHandler*)next())) {
298 handler->BeginEvent(entry);
299 }
300 next.Reset();
301 while((handler = (AliMCEventHandler*)next())) {
302 fMCEvent->AddSubsidiaryEvent(handler->MCEvent());
303 }
304 fMCEvent->InitEvent();
305 }
35152e4b 306
307 if (fPreReadMode == kLmPreRead) {
308 fMCEvent->PreReadAll();
309 }
310
93836e1b 311 return result;
312
5fe09262 313}
314
da97a08a 315void AliMCEventHandler::SelectParticle(Int_t i){
316 // taking the absolute values here, need to take care
317 // of negative daughter and mother
318 // IDs when setting!
93836e1b 319 if (TMath::Abs(i) >= AliMCEvent::BgLabelOffset()) i = fMCEvent->BgLabelToIndex(TMath::Abs(i));
320 if(!IsParticleSelected(TMath::Abs(i)))fParticleSelected.Add(TMath::Abs(i),1);
da97a08a 321}
322
323Bool_t AliMCEventHandler::IsParticleSelected(Int_t i) {
324 // taking the absolute values here, need to take
325 // care with negative daughter and mother
326 // IDs when setting!
327 return (fParticleSelected.GetValue(TMath::Abs(i))==1);
328}
329
330
331void AliMCEventHandler::CreateLabelMap(){
332
333 //
334 // this should be called once all selections where done
335 //
336
337 fLabelMap.Delete();
338 if(!fMCEvent){
339 fParticleSelected.Delete();
340 return;
341 }
342
343 VerifySelectedParticles();
da97a08a 344
345 Int_t iNew = 0;
93836e1b 346 for(int i = 0;i < fMCEvent->GetNumberOfTracks();++i){
da97a08a 347 if(IsParticleSelected(i)){
348 fLabelMap.Add(i,iNew);
349 iNew++;
350 }
351 }
352}
353
354Int_t AliMCEventHandler::GetNewLabel(Int_t i) {
93836e1b 355 // Gets the label from the new created Map
da97a08a 356 // Call CreatLabelMap before
357 // otherwise only 0 returned
358 return fLabelMap.GetValue(TMath::Abs(i));
359}
360
361void AliMCEventHandler::VerifySelectedParticles(){
362
363 //
364 // Make sure that each particle has at least it's predecessors
365 // selected so that we have the complete ancestry tree
366 // Private, should be only called by CreateLabelMap
367
368 if(!fMCEvent){
93836e1b 369 fParticleSelected.Delete();
370 return;
da97a08a 371 }
da97a08a 372
93836e1b 373 Int_t nprim = fMCEvent->GetNumberOfPrimaries();
374
375 for(int i = 0;i < fMCEvent->GetNumberOfTracks(); ++i){
376 if(i < nprim){
377 SelectParticle(i);// take all primaries
378 continue;
379 }
380
381 if(!IsParticleSelected(i))continue;
382
63a1afff 383 AliMCParticle* mcpart = (AliMCParticle*) fMCEvent->GetTrack(i);
93836e1b 384 Int_t imo = mcpart->GetMother();
385 while((imo >= nprim)&&!IsParticleSelected(imo)){
386 // Mother not yet selected
387 SelectParticle(imo);
0f0c06de 388 mcpart = (AliMCParticle*) fMCEvent->GetTrack(imo);
93836e1b 389 imo = mcpart->GetMother();
390 }
391 // after last step we may have an unselected primary
da97a08a 392 // mother
393 if(imo>=0){
394 if(!IsParticleSelected(imo))
395 SelectParticle(imo);
396 }
397 }// loop over all tracks
398}
399
5fe09262 400Int_t AliMCEventHandler::GetParticleAndTR(Int_t i, TParticle*& particle, TClonesArray*& trefs)
401{
402 // Retrieve entry i
c8b7b5d3 403 if (!fInitOk) {
404 return 0;
405 } else {
406 return (fMCEvent->GetParticleAndTR(i, particle, trefs));
407 }
5fe09262 408}
409
2d8f26f6 410void AliMCEventHandler::DrawCheck(Int_t i, Int_t search)
5fe09262 411{
412 // Retrieve entry i and draw momentum vector and hits
415d9f5c 413 fMCEvent->DrawCheck(i, search);
5fe09262 414}
415
890126ab 416Bool_t AliMCEventHandler::Notify(const char *path)
5fe09262 417{
53faeca4 418 // Notify about directory change
419 // The directory is taken from the 'path' argument
420 // Reconnect trees
421 TString fileName(path);
422 if(fileName.Contains("AliESDs.root")){
423 fileName.ReplaceAll("AliESDs.root", "");
424 }
4a1bbbb7 425 else if(fileName.Contains("AliESDs_wSDD.root")){
426 fileName.ReplaceAll("AliESDs_wSDD.root", "");
427 }
c8b5ce3a 428 else if(fileName.Contains("AliAOD.root")){
429 fileName.ReplaceAll("AliAOD.root", "");
430 }
53faeca4 431 else if(fileName.Contains("galice.root")){
432 // for running with galice and kinematics alone...
433 fileName.ReplaceAll("galice.root", "");
434 }
39734d21 435 else if (fileName.BeginsWith("root:")) {
436 fileName.Append("?ZIP=");
437 }
93836e1b 438
53faeca4 439 *fPathName = fileName;
93836e1b 440 AliInfo(Form("Path: -%s-\n", fPathName->Data()));
53faeca4 441
f6065654 442 ResetIO();
443 InitIO("");
444
93836e1b 445// Handle subsidiary handlers
446 if (fSubsidiaryHandlers) {
447 TIter next(fSubsidiaryHandlers);
448 AliMCEventHandler *handler;
449 while((handler = (AliMCEventHandler*) next())) {
450 TString* spath = handler->GetInputPath();
451 if (spath->Contains("merged")) {
452 if (! fPathName->IsNull()) {
453 handler->Notify(Form("%s/../.", fPathName->Data()));
454 } else {
455 handler->Notify("../");
456 }
457 }
458 }
459 }
460
5fe09262 461 return kTRUE;
462}
36e82a52 463
5fe09262 464void AliMCEventHandler::ResetIO()
465{
5efedd31 466// Clear header and stack
b72029a3 467
aff49450 468 if (fInitOk) fMCEvent->Clean();
5efedd31 469
36e82a52 470// Delete Tree E
415d9f5c 471 delete fTreeE; fTreeE = 0;
b72029a3 472
5efedd31 473// Reset files
89f249fc 474 if (fFileE) {delete fFileE; fFileE = 0;}
475 if (fFileK) {delete fFileK; fFileK = 0;}
476 if (fFileTR) {delete fFileTR; fFileTR = 0;}
19cfde04 477 fExtension="";
c8b7b5d3 478 fInitOk = kFALSE;
93836e1b 479
480 if (fSubsidiaryHandlers) {
481 TIter next(fSubsidiaryHandlers);
482 AliMCEventHandler *handler;
483 while((handler = (AliMCEventHandler*)next())) {
484 handler->ResetIO();
485 }
486 }
487
5fe09262 488}
489
490
491Bool_t AliMCEventHandler::FinishEvent()
492{
5efedd31 493 // Clean-up after each event
494 delete fDirTR; fDirTR = 0;
495 delete fDirK; fDirK = 0;
aff49450 496 if (fInitOk) fMCEvent->FinishEvent();
93836e1b 497
498 if (fSubsidiaryHandlers) {
499 TIter next(fSubsidiaryHandlers);
500 AliMCEventHandler *handler;
501 while((handler = (AliMCEventHandler*)next())) {
502 handler->FinishEvent();
503 }
504 }
505
5fe09262 506 return kTRUE;
507}
508
509Bool_t AliMCEventHandler::Terminate()
510{
511 // Dummy
512 return kTRUE;
513}
514
515Bool_t AliMCEventHandler::TerminateIO()
516{
517 // Dummy
518 return kTRUE;
519}
520
0a05cd41 521
0931e76a 522void AliMCEventHandler::SetInputPath(const char* fname)
0a05cd41 523{
524 // Set the input path name
525 delete fPathName;
526 fPathName = new TString(fname);
527}
93836e1b 528
529void AliMCEventHandler::AddSubsidiaryHandler(AliMCEventHandler* handler)
530{
531 // Add a subsidiary handler. For example for background events
532
533 if (!fSubsidiaryHandlers) fSubsidiaryHandlers = new TList();
534 fSubsidiaryHandlers->Add(handler);
535}