Warnings removed (R.Shahoyan)
[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
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];
195 sprintf(folder, "Event%d", iev);
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{
da97a08a 270 fParticleSelected.Delete();
271 fLabelMap.Delete();
9aea8469 272 // Read the next event
93836e1b 273
274 if (fEventsInContainer != 0) {
11e4d716 275 entry = (Long64_t) ( entry * Float_t(fNEvent) / Float_t (fEventsInContainer));
93836e1b 276 }
277
278
ed97dc98 279 if (entry == -1) {
280 fEvent++;
281 entry = fEvent;
282 } else {
283 fEvent = entry;
284 }
285
286 if (entry >= fNEvent) {
73bbf779 287 AliWarning(Form("AliMCEventHandler: Event number out of range %5lld %5d\n", entry, fNEvent));
9aea8469 288 return kFALSE;
289 }
93836e1b 290
291 Bool_t result = GetEvent(entry);
292
293 if (fSubsidiaryHandlers) {
294 TIter next(fSubsidiaryHandlers);
295 AliMCEventHandler *handler;
296 while((handler = (AliMCEventHandler*)next())) {
297 handler->BeginEvent(entry);
298 }
299 next.Reset();
300 while((handler = (AliMCEventHandler*)next())) {
301 fMCEvent->AddSubsidiaryEvent(handler->MCEvent());
302 }
303 fMCEvent->InitEvent();
304 }
35152e4b 305
306 if (fPreReadMode == kLmPreRead) {
307 fMCEvent->PreReadAll();
308 }
309
93836e1b 310 return result;
311
5fe09262 312}
313
da97a08a 314void AliMCEventHandler::SelectParticle(Int_t i){
315 // taking the absolute values here, need to take care
316 // of negative daughter and mother
317 // IDs when setting!
93836e1b 318 if (TMath::Abs(i) >= AliMCEvent::BgLabelOffset()) i = fMCEvent->BgLabelToIndex(TMath::Abs(i));
319 if(!IsParticleSelected(TMath::Abs(i)))fParticleSelected.Add(TMath::Abs(i),1);
da97a08a 320}
321
322Bool_t AliMCEventHandler::IsParticleSelected(Int_t i) {
323 // taking the absolute values here, need to take
324 // care with negative daughter and mother
325 // IDs when setting!
326 return (fParticleSelected.GetValue(TMath::Abs(i))==1);
327}
328
329
330void AliMCEventHandler::CreateLabelMap(){
331
332 //
333 // this should be called once all selections where done
334 //
335
336 fLabelMap.Delete();
337 if(!fMCEvent){
338 fParticleSelected.Delete();
339 return;
340 }
341
342 VerifySelectedParticles();
da97a08a 343
344 Int_t iNew = 0;
93836e1b 345 for(int i = 0;i < fMCEvent->GetNumberOfTracks();++i){
da97a08a 346 if(IsParticleSelected(i)){
347 fLabelMap.Add(i,iNew);
348 iNew++;
349 }
350 }
351}
352
353Int_t AliMCEventHandler::GetNewLabel(Int_t i) {
93836e1b 354 // Gets the label from the new created Map
da97a08a 355 // Call CreatLabelMap before
356 // otherwise only 0 returned
357 return fLabelMap.GetValue(TMath::Abs(i));
358}
359
360void AliMCEventHandler::VerifySelectedParticles(){
361
362 //
363 // Make sure that each particle has at least it's predecessors
364 // selected so that we have the complete ancestry tree
365 // Private, should be only called by CreateLabelMap
366
367 if(!fMCEvent){
93836e1b 368 fParticleSelected.Delete();
369 return;
da97a08a 370 }
da97a08a 371
93836e1b 372 Int_t nprim = fMCEvent->GetNumberOfPrimaries();
373
374 for(int i = 0;i < fMCEvent->GetNumberOfTracks(); ++i){
375 if(i < nprim){
376 SelectParticle(i);// take all primaries
377 continue;
378 }
379
380 if(!IsParticleSelected(i))continue;
381
63a1afff 382 AliMCParticle* mcpart = (AliMCParticle*) fMCEvent->GetTrack(i);
93836e1b 383 Int_t imo = mcpart->GetMother();
384 while((imo >= nprim)&&!IsParticleSelected(imo)){
385 // Mother not yet selected
386 SelectParticle(imo);
0f0c06de 387 mcpart = (AliMCParticle*) fMCEvent->GetTrack(imo);
93836e1b 388 imo = mcpart->GetMother();
389 }
390 // after last step we may have an unselected primary
da97a08a 391 // mother
392 if(imo>=0){
393 if(!IsParticleSelected(imo))
394 SelectParticle(imo);
395 }
396 }// loop over all tracks
397}
398
5fe09262 399Int_t AliMCEventHandler::GetParticleAndTR(Int_t i, TParticle*& particle, TClonesArray*& trefs)
400{
401 // Retrieve entry i
c8b7b5d3 402 if (!fInitOk) {
403 return 0;
404 } else {
405 return (fMCEvent->GetParticleAndTR(i, particle, trefs));
406 }
5fe09262 407}
408
2d8f26f6 409void AliMCEventHandler::DrawCheck(Int_t i, Int_t search)
5fe09262 410{
411 // Retrieve entry i and draw momentum vector and hits
415d9f5c 412 fMCEvent->DrawCheck(i, search);
5fe09262 413}
414
890126ab 415Bool_t AliMCEventHandler::Notify(const char *path)
5fe09262 416{
53faeca4 417 // Notify about directory change
418 // The directory is taken from the 'path' argument
419 // Reconnect trees
420 TString fileName(path);
421 if(fileName.Contains("AliESDs.root")){
422 fileName.ReplaceAll("AliESDs.root", "");
423 }
c8b5ce3a 424 else if(fileName.Contains("AliAOD.root")){
425 fileName.ReplaceAll("AliAOD.root", "");
426 }
53faeca4 427 else if(fileName.Contains("galice.root")){
428 // for running with galice and kinematics alone...
429 fileName.ReplaceAll("galice.root", "");
430 }
39734d21 431 else if (fileName.BeginsWith("root:")) {
432 fileName.Append("?ZIP=");
433 }
93836e1b 434
53faeca4 435 *fPathName = fileName;
93836e1b 436 AliInfo(Form("Path: -%s-\n", fPathName->Data()));
53faeca4 437
f6065654 438 ResetIO();
439 InitIO("");
440
93836e1b 441// Handle subsidiary handlers
442 if (fSubsidiaryHandlers) {
443 TIter next(fSubsidiaryHandlers);
444 AliMCEventHandler *handler;
445 while((handler = (AliMCEventHandler*) next())) {
446 TString* spath = handler->GetInputPath();
447 if (spath->Contains("merged")) {
448 if (! fPathName->IsNull()) {
449 handler->Notify(Form("%s/../.", fPathName->Data()));
450 } else {
451 handler->Notify("../");
452 }
453 }
454 }
455 }
456
5fe09262 457 return kTRUE;
458}
36e82a52 459
5fe09262 460void AliMCEventHandler::ResetIO()
461{
5efedd31 462// Clear header and stack
b72029a3 463
aff49450 464 if (fInitOk) fMCEvent->Clean();
5efedd31 465
36e82a52 466// Delete Tree E
415d9f5c 467 delete fTreeE; fTreeE = 0;
b72029a3 468
5efedd31 469// Reset files
89f249fc 470 if (fFileE) {delete fFileE; fFileE = 0;}
471 if (fFileK) {delete fFileK; fFileK = 0;}
472 if (fFileTR) {delete fFileTR; fFileTR = 0;}
19cfde04 473 fExtension="";
c8b7b5d3 474 fInitOk = kFALSE;
93836e1b 475
476 if (fSubsidiaryHandlers) {
477 TIter next(fSubsidiaryHandlers);
478 AliMCEventHandler *handler;
479 while((handler = (AliMCEventHandler*)next())) {
480 handler->ResetIO();
481 }
482 }
483
5fe09262 484}
485
486
487Bool_t AliMCEventHandler::FinishEvent()
488{
5efedd31 489 // Clean-up after each event
490 delete fDirTR; fDirTR = 0;
491 delete fDirK; fDirK = 0;
aff49450 492 if (fInitOk) fMCEvent->FinishEvent();
93836e1b 493
494 if (fSubsidiaryHandlers) {
495 TIter next(fSubsidiaryHandlers);
496 AliMCEventHandler *handler;
497 while((handler = (AliMCEventHandler*)next())) {
498 handler->FinishEvent();
499 }
500 }
501
5fe09262 502 return kTRUE;
503}
504
505Bool_t AliMCEventHandler::Terminate()
506{
507 // Dummy
508 return kTRUE;
509}
510
511Bool_t AliMCEventHandler::TerminateIO()
512{
513 // Dummy
514 return kTRUE;
515}
516
0a05cd41 517
0931e76a 518void AliMCEventHandler::SetInputPath(const char* fname)
0a05cd41 519{
520 // Set the input path name
521 delete fPathName;
522 fPathName = new TString(fname);
523}
93836e1b 524
525void AliMCEventHandler::AddSubsidiaryHandler(AliMCEventHandler* handler)
526{
527 // Add a subsidiary handler. For example for background events
528
529 if (!fSubsidiaryHandlers) fSubsidiaryHandlers = new TList();
530 fSubsidiaryHandlers->Add(handler);
531}